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A transformation  technique  is  introduced,  by  means  of 
which  the  inequality  constrained  problem  is  converted  into  an 
equality  constrained  problem  involving  differential  constraints, 
terminal  constraints,  and  a control  equality  constraint.  The 
transformation  technique  takes  advantage  of  the  partial  linearity  of 
the  state  inequality  constraint  so  as  to  yield  a transformed 
problem  characterized  by  a new  state  vector  of  minimal  size. 

This  concept  is  important  computationally,  in  that  the  com- 
puter time  per  iteration  increases  with  the  square  of  the 
dimension  of  the  state  vector. 

In  order  to  illustrate  the  advantages  of  the  new  trans- 
formation technique,  several  numerical  examples  are  solved  by 
means  of  the  sequential  gradient-restoration  algorithm  for 
optimal  control  problems  involving  nondifferential  constraints. 
The  examples  show  the  substantial  savings  in  computer  time 
for  convergence,  which  are  associated  with  the  new  transforma- 
tion technique. 

Key  Words . Optimal  control,  numerical  methods,  computing 
methods,  transformation  techniques,  sequential  gradient- 
restoration  algorithm,  nondifferential  constraints,  state 
inequality  constraints,  linear  state  inequality  constraints, 
partially  linear  state  inequality  constraints. 
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1 . Introduction 

In  recent  years,  considerable  attention  has  been  devoted 
to  the  study  of  optimal  control  problems  with  bounded  state. 

The  approaches  employed  are  of  two  types:  (i)  the  direct  ap- 
proach, in  which  a predetermined  number  and  sequence  of  subarcs 
are  assumed,  and  (ii)  the  indirect  approach  in  which  no  pre- 
determined number  and  sequence  of  subarcs  are  assumed. 

Direct  Approach.  In  the  direct  approach,  the  extremal  arc 
is  viewed  as  being  composed  of  several  subarcs,  some  internal 
to  the  state  boundary  and  some  lying  on  the  state  boundary.  At 
the  points  of  junction  of  different  subarcs,  discontinuities  in 
the  control  and  the  multipliers  are  allowed.  In  practice,  this 
view  of  problems  with  state  inequality  constraints  can  be  im- 
plemented through  a variety  of  techniques,  all  of  them  requir- 
ing satisfaction  of  the  state  inequality  constraint,  with  this 
understanding:  strict  inequality  must  be  enforced  for  the 
subarcs  internal  to  the  state  boundary,  and  strict  equality 
must  be  enforced  for  the  subarcs  lying  on  the  state  boundary. 

Indirect  Approach.  In  the  indirect  approach,  the  extremal 
arc  is  viewed  as  a single  subarc,  even  though  a portion  of  it 
may  lie  extremely  close  to  the  state  boundary.  Along  this 
single  subarc,  the  control  and  the  multipliers  are  regarded  to 
be  continuous  functions  of  the  time.  In  practice,  this  view  of 
problems  with  state  inequality  constraints  can  be  implemented 
through  either  penalty  function  techniques  or  transformat  ion 
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techniques.  The  latter  techniques  (see,  for  instance,  Refs. 
1-5)  constitute  the  key  concern  of  this  paper. 

Transformation  Techniques.  In  a transformation  scheme,  a 
Valentine- type  representation  is  employed  for  the  state  ine- 
quality constraint  Lq(x,  0 ) >0  (Ref.  4).  If  k is  the  order  of 

5 

the  state  inequality  constraint  , the  dimension  of  the  state 
vector  is  increased  by  k.  That  is,  the  n-vector  state  x is 
augmented  by  the  k- vector  y,  which  is  defined  in  such  a way 
that  the  state  inequality  constraint  L^(x, 0 ) > 0 is  satisfied 
automatically  at  each  point  of  the  trajectory. 

Concerning  the  m-vector  control  u,  one  has  two  options: 

(i)  to  leave  its  dimension  unchanged  through  a substitution 
technique  (Ref.  4);  or  (ii)  to  increase  its  dimension  by  one. 
In  the  second  option,  the  m-vector  control  u is  augmented  by 
the  scalar  control  w.  In  turn,  the  augmented  control  (u,w) 
must  satisfy  a nondifferential  constraint  everywhere  along  the 
trajectory  (Ref.  5) . 

Research  Objectives.  A drawback  of  the  transformation 
schemes  of  Refs.  4 and  5 is  that  they  cause  an  increase  in  the 


A state  inequality  constraint  Lq(x , 0 ) £ 0 is  defined  to  be  of 
order  k if  the  kth  time  derivative  of  the  function  Lq(x,6) 
is  the  first  to  contain  the  control  explicitly. 
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dimension  of  the  state  vector  and,  consequently,  an  increase 
in  the  computer  time  per  iteration.  As  an  example,  if  one 
uses  the  sequential  gradient-restoration  algorithm  for  optimal 
control  problems  with  nondifferential  constraints  (Ref.  5), 
the  computer  time  per  iteration  is  proportioned  to  the  square 
of  the  dimension  of  the  state  vector.  Clearly,  in  devising 
transformation  techniques  converting  an  inequality  constrained 
problem  into  an  equivalent  equality  constrained  problem,  it  is 
important  to  prevent  or  limit  the  increase  in  the  dimension  of 
the  state  vector;  in  other  words,  it  is  important  that  the 
dimension  of  the  new  state  vector  be  as  small  as  possible. 

In  this  paper,  we  pursue  this  point  of  view  with  reference 
to  a special  category  of  optimal  control  problems,  namely, 
problems  where  the  state  inequality  constraint  Lq(x,0)>O  is 
linear  in  some  or  all  of  the  components  of  the  state  vector. 

We  consider  the  k-vector  function  L(x,0),  composed  of  Lq(x,0) 
and  its  first  k-1  derivatives.  We  make  the  following  assump- 
tions: (i)  the  n-vector  state  x can  be  partitioned  into  vec- 

tors x^  and  xB  having  dimensions  a and  b,  respectively;  (ii) 

the  k-vector  function  L(xa,xd,0)  can  be  partitioned  into  two 

A D 

vector  functions  Kx,  ,xn,8)  and  L(xa,x_,0)  having  dimensions 
a and  c,  respectively;  and  (iii)  the  a-vector  function  L(xa,x0,0) 
is  linear  in  x . With  this  understanding,  the  elimination  of 

A 

x^  becomes  possible  and  one  can  formulate  a new  optimal  control 


problem  involving  a state  vector  of  reduced  size.  In  other 
words,  the  augmented  state  vector  (x,y)  having  dimension  n + k 
is  replaced  by  the  new  state  vector  (xD,y)  having  dimension 
b + k,  with  the  following  implication:  if  the  algorithm  of 
Ref.  5 is  employed,  the  computer  time  per  iteration  is  scaled 
according  to  the  square  of  the  ratio  (b  + k)/(n  + k). 

In  theory,  the  transformation  technique  proposed  here 
appears  to  have  a very  particular  nature.  In  practice,  a large 
number  of  technical  problems  fall  within  the  scheme  described 
here.  With  this  in  mind,  the  advantages  of  the  new  transforma- 
tion technique  are  illustrated  through  several  examples,  which 
are  solved  by  means  of  the  sequential  gradient-restoration 
algorithm  for  optimal  control  problems  with  nondifferential 
constraints  (Ref.  5) . The  examples  show  the  substantial  savings 
in  computer  time  for  convergence,  which  are  associated  with 
the  new  transformation  technique. 
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2.  Statement  of  the  Problem 

2.1.  Notation.  Let  0 denote  the  independent  variable, 

and  let  x(6),  u(6)  denote  the  dependent  variables.  The  time  6 

is  a scalar,  the  state  x(0)  is  an  n-vector , and  the  control 

u ( 0 ) is  an  m-vector.6  The  initial  time  is  0 = 0 and  the  final 

7 

time  is  0 = x.  The  final  time  is  either  given  or  tree. 

2.2.  Optimization  Problem.  With  the  above  definitions, 
the  optimization  problem  can  be  stated  as  follows. 

Problem  Pi.  Minimize  the  functional 

1=  l f (x,  u,  0 ) d0  + [g  (x,  0 )]  T (1) 

Jo 

with  respect  to  the  state  x(6),  the  control  u(6),  and  the  para- 
meter x which  satisfy  the  differential  constraints 

dx/d©  = 4>  (x,u,  0) , 0 < 0 < x , (2) 

the  state  inequality  constraint 

Lq(x,  0 ) >0,  0 < 0 < x , (3) 

and  the  boundary  conditions 


6A11  vectors  are  column  vectors. 

7In  the  latter  case,  x is  a parameter  to  be  optimized. 


0(x,e)]T  = o . (5) 


In  the  above  equations,  the  quantities  I , f , g , Lq  are  scalar , the 
function  <J>  is  an  n-vector,  and  the  function  ip  is  a q-vector. 
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3 . Transformation  of  the  Problem 

In  this  section,  we  convert  the  inequality  constrained 
problem  into  an  equality  constrained  problem.  We  employ  the 
Valentine- type  transf ormation  proposed  by  Jacobson  in  Ref.  4. 

We  assume  that  the  state  inequality  constraint  (3)  has  order  k. 
This  means  that  the  kth  time  derivative  of  the  function  Lg(x , 0 ) 
is  the  first  to  contain  the  control  u explicitly. 

3.1.  Valentine-Type  Transformation.  We  introduce  the 
auxiliary  state  vector  y(G)  and  the  auxiliary  control  variable 
w(0),  where  y is  a k-vector  and  w is  a scalar.  We  require  the 
k-vector  y and  the  scalar  w to  satisfy  the  differential  equa- 
tions 

dv1/d0  = y2  , 
dy2/d0  = y3  , 

(6) 

dyk_1/d0  = yk, 

dyk/dO  = w . 

Next,  we  discard  the  state  inequality  constraint  (3)  and  re- 
place it  with  the  state  equality  constraint 


Lq  (x,0)  - Aq  (y)  =9, 


a0  (y)=  * 


(7) 
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where  y^(0)  denotes  the  first  component  of  the  auxiliary  state 
vector  y(9).  Then,  we  execute  k time  derivatives  of  both 
sides  of  Eq.  (7).  This  yields  the  relations 


L1  (x , 0 ) - A^  (y)  = 0 , 
l2 (x, 6)  - a2 (y)  = 0 , 
L3 (x, 6 ) - A3 (y)  =0  , 
L4  ( x » 0 ) - A4 (y)  = 0 , 


Lk-1 (X' 6)  “ Ak-1 (y)  = 0 ' 


Mk(x,u,0)  - Bk(y,w)  = 0.  (9) 

In  Eqs . (8) - (9) , the  symbols  1^  and  denote  the  successive  total 

derivatives  of  Lq(x,0)  , and  the  symbols  A^  and  denote  the  successive 
total  derivatives  of  AQ(y).  Specifically,  the  derivatives 
A^  (y) , i = 0,1,... ,k-l,  are  given  by 


t 


A0 (y)  = Y1  ' 

(y)  = 2yxy2  , 

A2(y)  =2y2  + 2Yly3  , 

A3 (y)  = 6Y2^3  + 2^1^4  ' 

A4  (y)  = 6y3  + 8y2y4  + 2y;Ly5  , 


(10) 


(ID 


and  the  derivative  B (y,w)  is  given  by 

K 


Bx (y ,w)  = 2y3w, 

B2 (y , w)  = 2y2  + 2yxw  , 

B3 (y , w)  = 6y2y3  + 2y3w  , 

B4  (y,w)  = 6y3  + 8y2y4  + 2yxw  , if  k=4, 

B5(y,w)  = 20y3y4  + 10y2y5  + 2y3w, 

Note  that  Eqs.  (7)  and  (8)  are  successive  first  integrals 
of  Eq . (9).  Also  note  that  Eqs.  (7) -(8),  when  applied  at  the 

initial  point,  yield  the  initial  conditions  for  the  auxiliary 
state  vector.  Since  x(0)  is  given,  Eqs.  (7) -(8)  imply  that 


if 

k=  1, 

if 

k = 2, 

if 

k = 3, 

if 

ii 

if  k = 

y ( 0 ) = given 


(12) 


L 
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Summarizing,  the  state  inequality  constraint  (3)  can  be 
replaced  by  the  state  equality  constraint  (7).  In  turn,  (7) 
can  be  replaced  by  the  control  equality  constraint  (9),  which 
is  to  be  employed  in  combination  with  the  differential  equa- 
tions (6)  and  the  initial  conditions  (12) . 

3.2.  Notation.  In  more  compact  notation,  the  differen- 
tial equations  (6)  can  be  rewritten  in  the  form 


dy/dB  = to(y , w)  , 0 < 9<  x , 

where  y and  denote  the  k-vectors 

y = [yi'y2 yk-i'yklT' 

T 

co(y,w)  = [y2,y3, »Yk»w]  • 

In  turn,  the  relations  (7)- (8)  can  be  rewritten  as 

L(x,0)  - A (y ) = 0 , 

where  L(x,6)  and  A(y)  denote  the  k-vectors 

L (x,  0 ) = [Lq(x,0)  ,L^(x,0)  ,Lk-l  (X/  e )]  T ' 

A(y)  = [AQ  (y)  ,A1  (y)  'Ak-1  (y)JT  • 

Finally,  the  relation  (9)  can  be  rewritten  as 

S (x,y,u,w,Q)  = 0 , 


(13) 


(14) 


(15! 


(16) 


(17: 


where 


S(x,y,u,w,0)  =Mk(x,u,6)  -Bk(y,w)  . (18) 

With  the  aid  of  this  notation,  Problem  PI  can  be  reformulated 
as  follows. 

Problem  P2.  Minimize  the  functional 

r 

1=  l f (x,u,  0)d6  + [ g (x,  0)]  (19) 

J 0 

with  respect  to  the  state  x(0),  y(0),  the  control  u(0),w(0), 
and  the  parameter  x which  satisfy  the  differential  constraints 


dx/d0  = <f>  (x,u,  0 ) , 

0 < 0 < T , 

(20) 

dy/d0  = m (y,w) , 

0 < 0 < T , 

(21) 

the 

nondifferential  constraint 

S (x,y,u,w,0)  = 0, 

0 < 0 < T , 

(22) 

and 

the  boundary  conditions 

x ( 0 ) = given , 

y (0)  = given. 

(23) 

[>  (x,0)]T 

= 0 . 

(24) 

In  the  above  equations,  the  quantities  I,f,g,S  are  scalar,  the 
function  is  an  n-vector,  the  function  ■ is  a k-vector,  and 
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the  function  ip  is  a q-vector. 

Note  that  the  vectors  x(6)  and  y(0)  are  automatically 
consistent  with  the  relation 

L(x,0)  - A (y ) =0,  0 < 6 < x,  (25) 

at  every  point  of  the  trajectory,  providing  they  are  consistent 
with  (25)  at  the  initial  point.  Therefore,  the  initial  condi- 
tions (23)  must  be  such  that 

L (x (0) ,0)  - A(y (0) ) = 0 . (26) 

3.3.  Partitioning  of  the  State  Vector.  Assume  now  that  the 

n-vector  state  x(0)  is  partitioned  into  vectors  x (0)  and  x (f* ) 

A B 

having  dimensions  a and  b,  respectively,  such  that 

a + b = n,  0<a<k,  n-k  < b < n . (27) 

Under  these  conditions.  Problem  P2  can  be  reformulated  as  fol- 
lows . 

Problem  P2.  Minimize  the  functional 

f 

1=  \ f (xA,xB,u,0)dO  + [g(xA,xB,0)]  t (28) 

Jo 

with  respect  to  the  state  x.(G),  xn(0),  y(0),  the  control  u(n), 

A B 

w(Q),  and  the  parameter  x which  satisfy  the  differential  con- 


straints 
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A A A B 

dVd0  = VXA'XB' 

dy/d0  = w (y , w) , 

the  nondifferential  constraint 

S(XA'XB'Y'U'W‘ 
and  the  boundary  conditions 


xA(0)  = given, 


B 


consistent  with  the  relation 


A'  B' 


u,  6)  , 

0 < 0 

1 T ' 

(29) 

u,0)  , 

0 < 0 

5 1 » 

(30) 

CD 

V 1 

o 

< T ' 

(31) 

= o, 

0 < 0 < 

T , 

(32) 

= given , 

y(o: 

) = given , 

(33) 

• 

(34) 

antities 

iff 

3 are 

i scalar, 

the 

: function 

4>b  is 

a b- 

■vector , 

the 

the  funct 

ion  iJj  : 

Ls  a 

q-vector 

• 

6) , xB(0) 

, y(0) 

are 

automatically 

= o , 

0 < 0 < 

T , 

(35) 

at  any  point  of  the  trajectory,  providing  they  are  consistent 
with  (35)  at  the  initial  point.  Therefore,  the  initial  con- 
ditions (33)  must  be  such  that 


L (xA(0)  ,xR(0)  ,0)  - A ( y (0))  = 0. 


'B 


(36) 


J 
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Full  advantage  of  the  partitioning  of  the  state  vector  x 

into  vectors  x-  and  xD  is  taken  in  the  following  section. 

A fcj 
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4 . Further  Transformation 

In  many  technical  applications,  the  left-hand  side  of  the 
state  inequality  constraint  Lq(x,9)  is  linear  in  some  or  all 
of  the  components  of  the  state  vector.  In  addition,  depending 
on  the  form  of  the  differential  system,  some  of  the  successive 
total  derivatives  L^(x,0)  might  be  linear  in  some  or  all  of  the 
components  of  the  state  vector.  In  other  words,  it  is  plausible 
that  situations  exist  such  that  the  following  properties  hold. 

(i)  The  n-vector  state  x(0)  can  be  partitioned  into  vec- 
tors XA(0)  and  xg(9)  having  dimensions  a and  b,  respectively, 
such  that 

a + b = n,  0 < a < k,  n - k < b < n . (37) 

(ii)  The  k-vector  function  L(x:.,xc>,0)  can  be  partitioned 

A 15 

~ A 

into  two  vector  functions  L(xA,xB,0)  and  L(xA,Xg,0)  having 
dimensions  a and  c,  respectively,  such  that 

a + c = k , 0 < a < k , 0<c<k.  (38) 

(iii)  The  a-vector  function  L(xA,Xg,0)  is  linear  in  x^, 
that  is, 

L (xA,xB, 0 ) = M (0) xA + N (xB, 0 ) , (39) 

where  M(9)  is  an  ax  a matrix  and  N(x  ,0)  is  an  a-vector. 

Under  the  above  conditions,  important  simplifications  are 
possible.  If  the  square  matrix  M(0)  is  nonsingular,  Eq . (39) 


I 
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admits  the  formal  solution 

xA=  M_1  (0)[L(xA,xB,0)- N(xb,6)]  , (40) 

which  in  the  light  of  (35)  can  be  rewritten  as 

xA=  M_1(0)[A(y)- N(xb,0)]  . (41) 

Here,  A(y)  is  obtained  by  partitioning  the  k-vector  function 
A(y)  into  two  vector  functions  A(y)  and  A(y)  having  dimensions 
a and  c,  respectively. 

Under  the  above  premises,  the  elimination  of  xA  becomes 
possible,  and  one  can  formulate  a new  optimal  control  problem 
involving  the  reduced  state  vector  xn  as  well  as  the  auxiliary 
state  vector  y. 

Problem  P3.  Minimize  the  functional 

1=1  f (xB,y,u,0)dO  + [g(xB,y,e)]T  (42) 

J 0 

with  respect  to  the  state  x (6),y(0),  the  control  u(0),w(0), 
and  the  parameter  x which  satisfy  the  differential  constraints 


dxB/de  = <|.B(xB,y,u,e), 

0 < 0 < x , 

(43) 

dy/d6  = a)  (y  ,w)  , 

0 < 0 < x , 

(44) 

the  nondifferential  constraint 

S(x0,y,u,w,O)  = 0,  0 < 0 < t , 


(45) 
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and  the  boundary  conditions 

xD (0)  = given,  y (0)  = given,  (46) 

D 

[iMxB,y,0)]T  = 0 . (47) 

In  the  above  equations,  the  quantities  I,f,g,S  are  scalar,  the 

function  <J>_  is  a b-vector,  the  function  w is  a k-vector,  and 
B 

g 

the  function  iJj  is  a q-vector. 

4.1.  Remark.  Comparison  of  Problem  P2  with  Prob- 
lem P3  shows  that  the  augmented  state  vector  (x,y)  having  di- 
mension n + k has  been  replaced  by  the  new  state  vector  (x  ,y) 
having  dimension  b + k,  with  the  following  implication:  if  the 
algorithm  of  Ref.  5 is  employed,  the  computer  time  per  iteration 
is  scaled  according  to  the  square  of  the  ratio  (b  + k)/ (n  + k). 

In  theory,  the  transformation  technique  proposed  here  ap- 
pears to  have  a very  particular  nature.  In  practice,  a large 
number  of  technical  problems  fall  within  the  scheme  described 
here.  With  this  in  mind,  the  advantages  of  the  new  transforma- 
tion technique  are  illustrated  in  the  following  sections  through 

••  m 

several  examples,  which  are  solved  by  means  of  the  sequential 
gradient-restoration  algorithm  for  optimal  control  problems 
with  nondifferential  constraints  (Ref.  5).  The  examples  show 
the  substantial  savings  in  computer  time  for  convergence,  which 
are  associated  with  the  new  transformation  technique. 

O 

After  Problem  P3  is  solved  and  the  functions  xb(')»  y(0)  are 
known,  the  function  x ( ) can  be  computed  a posterior i with (41). 


— 
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5 . Description  of  the  Algorithm 

With  the  sequential  ordinary  gradient-restoration  algori- 
thm for  optimal  control  problems  involving  nondifferential 
constraints  (SOGRA,  Refs.  5-6),  numerical  solutions  can  be 
obtained  using  either  the  formulation  of  Problem  P2  or  the 
formulation  of  Problem  P3. 

5.1.  Time  Normalization.  The  sequential  ordinary  gradient- 
restoration  algorithm  (SOGRA)  requires  that  the  actual  time  6 
be  replaced  by  the  normalized  time 

t = 8/t  , 0 < t < 1 , (48) 

defined  in  such  a way  that  t = 0 at  the  initial  point  and  t=  1 

at  the  final  point.  The  actual  final  time  t,  if  it  is  free, 

9 

becomes  a parameter  to  be  optimized.  In  this  way,  an  optimal 
control  problem  with  variable  final  time  is  converted  into  an 
optimal  control  problem  with  fixed  final  time.  Even  though 
this  replacement  of  independent  variable  is  necessary  computa- 
tionally, it  was  not  employed  in  the  description  of  Problems 
PI  through  P3,  for  the  sake  of  brevity.  However,  this 


9 

The  dimension  of  the  parameter  is  p = 0 if  t is  fixed  and  p=  1 
if  t is  free. 
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replacement  of  independent  variable  is  employed  in  the 
description  of  the  numerical  examples  of  Section  7. 

5.2.  Auxiliary  Functionals.  In  addition  to  the  function- 
al I,  several  auxiliary  functionals  are  of  interest  in  the 
implementation  of  SOGRA.  These  auxiliary  functionals  are 
denoted  by  the  symbols  J,  P,  Q and  have  the  following  meaning: 

J is  the  augmented  functional,  that  is,  the  functional  I 
augmented  linearly  by  the  constraints  through  Lagrange  multi- 
pliers; P is  the  constraint  error,  that  is,  the  norm  squared 

of  the  error  in  the  feasibility  equations;  and  Q is  the  optim- 
ality condition  error,  that  is,  the  norm  squared  of  the  error 
in  the  optimality  conditions.  In  the  definitions  of  J,  P,  Q, 
it  is  tacitly  assumed  that  the  given  initial  conditions  are 
satisfied  at  every  iteration  of  SOGRA. 

5.3.  Convergence  Conditions.  The  auxiliary  functionals 
P and  Q are  defined  in  such  a way  that 

P=Q=0  (49) 

for  the  optimal  solution.  For  an  approximation  to  the  optimal 
solution,  Eqs.  (49)  are  replaced  by  the  inequalities 

P < , Q < c2  ' (50) 

where  and  are  small,  preselected  numbers. 


^ . 
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5.4.  Sequential  Ordinary  Gradient-Restoration  Algorithm. 
The  sequential  ordinary  gradient-restoration  algorithm  (SOGRA, 
Refs.  5-6)  is  an  iterative  technique  which  includes  a sequence 
of  cycles  having  the  following  properties:  (i)  the  functions 
available  both  at  the  beginning  and  at  the  end  of  each  cycle 
are  feasible;  that  is,  they  are  consistent  with  the  feasibility 
equations  within  the  preselected  accuracy  (50-1);  and  (ii)  the 
functions  produced  at  the  end  of  each  cycle  are  characterized 
by  a value  of  the  the  functional  I which  is  smaller  than  that 
associated  with  the  functions  available  at  the  beginning  of  the 
cycle. 

To  achieve  the  above  properties,  each  cycle  is  made  of 
two  phases,  a gradient  phase  and  a restoration  phase.  These 
phases  are  now  described. 

The  gradient  phase  is  started  only  when  the  constraint 
error  P satisfies  Ineq.  (50-1).  It  involves  a single  itera- 
tion, which  is  designed  to  decrease  the  value  of  the  functional 
I or  the  augmented  functional  J,  while  satisfying  the  con- 
straints to  first  order.  During  this  iteration,  the  first 
variation  of  the  functional  I is  minimized,  subject  to  the 
linearized  constraints  and  a quadratic  constraint  on  the  vari- 
ations of  the  control  and  the  parameter. 

The  restoration  phase  is  started  only  when  the  constraint 
error  P violates  Ineq.  (50-1).  The  restoration  phase  involves 


1 


one  or  more  iterations.  In  each  restorative  iteration,  the 
objective  is  to  reduce  the  constraint  error  P,  while  the  con- 
straints are  satisfied  to  first  order  and  the  norm  squared  of  the 
variations  of  the  control  and  the  parameter  is  minimized.  The  restor- 
ation phase  is  terminated  whenever  Ineq.  (50-1)  is  satisfied. 

In  summary,  the  main  properties  of  the  sequential  ordi- 
nary gradient-restoration  algorithm  can  be  written  as  follows: 


J < J , 

p < , 

(51) 

p < p , 

/\ 

P > , 

(52) 

I < I , 

p 1 ' 

P < £-1  , 

(53) 

where  (51)  hold  for  a gradient  iteration,  (52)  hold  for  a res- 
torative iteration,  and  (53)  hold  for  a complete  gradient- 
restoration  cycle.  In  the  above  relations,  I,  J,  P denote 
quantities  at  the  beginning  of  an  iteration  or  a cycle;  I,  j, 

P denote  quantities  at  the  end  of  an  iteration;  and  i,  J,  P 
denote  quantities  at  the  end  of  a cycle. 

Each  iteration,  gradient  or  restorative,  involves  two  dis- 
tinct operations:  (i)  the  computation  of  the  displacements  per 
unit  stepsize  and  (ii)  the  computation  of  the  stepsize  t.  Here, 
the  terminology  of  Refs.  5-6  is  employed.  Operation  (i)  entails 
the  solution  of  a linear,  two-point  boundary-value  problem  by 
means  of  the  method  of  particular  solutions  (Refs.  7-9).  Opera- 
tion (ii)  entails  the  execution  of  a one-dimensional  search 
(Refs.  10-11) . 
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5-5-  Computational  Effort.  At  each  iteration  of  the 


gradient  phase  or  the  restoration  phase,  a linear,  two-point 


boundary-value  problem  must  be  solved.  If  the  method  of  parti- 


cular solutions  is  employed  (Refs.  7-9) , one  must  execute 


n^+P+1  independent  sweeps  of  the  linearized  system  governing 


the  variations  associated  with  the  gradient  phase  or  the 


restoration  phase,  where  n^  is  the  dimension  of  the  state  vec- 


tor and  p is  the  dimension  of  the  parameter.  Note  that 


n2 = n + k for  Problem  P2, 


(54-1) 


n^ = b + k for  Problem  P3. 


(54-2) 


Hence,  the  following  factor  is  indicative  of  the  algorithmic 


time  per  iteration  : 


Wi  = 2ni  (ni  + p + 1), 


i = 2 or  3 , 


with  the  implication  that 


W2  = 2 (n  + k)  (n  + k + p + 1)  for  Problem  P2 


VJ  = 2 (b  + k)  (b  + k + p + 1)  for  Problem  P3. 


The  subscript  i = 2 is  employed  for  Problem  P2,  and  the  subs- 
cript i = 3 is  employed  for  Problem  P3. 
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5.6.  Gradient  Stepsize.  The  computation  of  the  gradient 
stepsize  requires  the  consideration  of  the  functions 

J = J (a)  , P = P (a)  . (58) 

Then,  a one-dimensional  search  scheme  is  applied  to  (58-1) , and 
a value  of  the  stepsize  a is  selected  for  which  the  following 
relations  are  satisfied: 

J (a)  < J (0) , P(a)<P*,  t (a)  > 0 , (59) 

where  x is  the  final  time  and  P*  is  a preselected  number,  not 
necessarily  small.  Satisfaction  of  Ineq.  (59-1)  is  possible 
because  of  the  descent  property  of  the  gradient  phase.  Ineq. 
(59-2)  is  introduced  to  prevent  excessive  constraint  violation. 
And  Ineq.  (59-3)  is  required  for  problems  with  free  final  time. 

Prior  to  the  satisfaction  of  (59) , a scanning  process  is 
employed,  leading  to  the  bracketing  of  the  minimum  point  for 
J(a).  This  operation  is  then  followed  by  a Hermitian  cubic 
interpolation  process  (Refs.  10-11) , which  is  stopped  whenever 
the  following  relation  is  satisfied:'*''*' 

|J  (a) | < e,  or  |j  (a)/J  (0) | < e , (60) 

1 a 1 - 3 1 ot  a - 4 


11 


rhe  symbols  ar>d 


4 


denote  small,  preselected  numbers. 
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subject  to  an  upper  limit  for  the  number  of  search  steps  N . 

Once  a stepsize  aQ  has  been  selected  consistently  with  either 
(60)  or  the  prescribed  upper  limit  for  the  number  of  search 
steps,  Ineqs. (59)  must  be  checked.  If  satisfaction  occurs, 
then  the  stepsize  aQ  is  accepted.  If  any  violation  occurs, 
then  the  stepsize  aQ  must  be  bisected  progressively  until  sat- 
isfaction of  (59)  is  finally  achieved. 

5.7.  Restoration  Stepsize.  The  computation  of  the  restor- 
ation stepsize  requires  the  consideration  of  the  function 

P = P(ct)  . (61) 

Then,  the  stepsize  a must  be  selected  so  that  the  following 
relations  are  satisfied: 

P (a)  < P (0)  , f (a)  > 0 . (62) 

Satisfaction  of  Ineq.  (62-1)  is  possible  because  of  the  des- 
cent property  of  the  restoration  phase.  Ineq.  (62-2)  is 
required  for  problems  with  free  final  time. 

In  order  to  achieve  satisfaction  of  (62) , a bisection 
process  is  applied  to  the  restoration  stepsi ze  < , starting  from 
the  reference  stepsize  a =1.  This  reference  stepsize  has  the 
property  of  yielding  one-step  restoration  for  the  special  case 
where  all  the  constraints  are  linear. 
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5.8.  Iterative  Procedure  for  the  Restoration  Phase.  The 
descent  property  (52-1)  of  the  restoration  phase  guarantees 
satisfaction  of  Ineq.  (62-1)  at  the  end  of  any  iteration,  but 
not  satisfaction  of  Ineq.  (50-1) . Therefore,  the  restoration 
algorithm  must  be  employed  iteratively  until  Ineq.  (50-1)  is 
satisfied.  At  this  point,  the  restoration  phase  is  terminated. 

5.9.  Descent  Property  of  a Cycle.  A descent  property 
exists  for  a complete  gradient-restoration  cycle  under  the 

~ /\ 

assumption  of  small  stepsizes.  More  specifically,  let  I,  I,  I 
denote  the  values  of  the  functional  under  consideration  at  the 
beginning  of  the  gradient  phase,  at  the  end  of  the  gradient 
phase,  and  at  the  end  of  the  subsequent  restoration  phase. 

Note  that  I and  I are  not  comparable,  since  the  constraints  are 
not  satisfied  to  the  same  accuracy.  On  the  other  hand,  I and  I 
are  comparable,  and  the  gradient  stepsize  can  be  selected  so 
that 

I < I . (63) 

This  inequality  constitutes  the  descent  property  of  a complete 
gradient-restoration  cycle.  In  order  to  enforce  it,  one  pro- 
ceeds as  follows.  At  the  end  of  the  restoration  phase,  one 
must  verify  Ineq.  (63).  If  it  is  satisfied,  the  next  gradient 
phase  is  started;  otherwise,  the  previous  gradient  stepsize  is 
bisected  as  many  times  as  needed  until,  after  restoration, 

Ineq.  (63)  is  satisfied. 


r 


26 


AAR-137 


6 . Experimental  Conditions 

In  order  to  evaluate  the  transformation  techniques  dis- 
cussed in  Sections  3-4,  several  numerical  examples  were  solved. 
The  sequential  gradient-restoration  algorithm  of  Ref.  5 was 
programmed  in  FORTRAN  IV,  and  the  numerical  results  were 
obtained  in  double-precision  arithmetic. 

Computations  were  performed  at  Rice  University  using  an 
IBM  370/155  computer.  For  each  example,  the  interval  of  in- 
tegration was  divided  into  50  steps.  The  differential  equations 
were  integrated  using  Hamming's  modified  predictor-corrector 
method  with  a special  Runge-Kutta  starting  procedure  (Ref.  12). 
The  definite  integrals  I,  J,  P,  Q were  computed  using  a modified 
Simpson's  rule.  The  method  of  particular  solutions  (Refs.  7-9) 
was  used  to  solve  the  linear,  two-point  boundary-value  problems 
associated  with  both  the  gradient  phase  and  the  restoration 
phase . 

6.1.  Convergence  Conditions.  The  parameters  e.,  e4 

appearing  in  Ineqs.  (50)  and  (60)  were  set  at  the  levels  -*-2 

e^E-08,  e2  = E-04,  e4  = E - 03  . (64) 

The  tolerance  level  (64-1)  characterizes  the  restoration  phase; 


19  4- 

The  symbol  E + ab  stands  for  10- 


t 
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the  tolerance  levels  (64-1)  and  (64-2)  , employed  in  combina- 
tion, characterize  the  algorithm  as  a whole;  and  the  tolerance 
level  (64-3)  characterizes  the  one-dimensional  search  for  the 
gradient  stepsize. 

6.2.  Safeguards . For  the  gradient  phase,  the  parameter 
P*  appearing  in  Ineq.  (59-2)  was  set  at  the  level 

P*  = 10  . (65) 


The  tolerance  level  (65)  limits  the  constraint  violation  which 
is  permissible  during  the  gradient  phase.  Also  for  the  grad- 
ient phase,  the  number  of  Hermitian  search  steps  Ng  required  to 
satisfy  Ineq.  (60)  was  subject  to  the  upper  bound 

N < 5 . (66) 

s - 

6.3.  Nonconvergence  Conditions.  The  sequential  gradient- 

restoration  algorithm  was  programmed  to  stop  whenever  violation 

1 3 

of  any  of  the  following  inequalities  occurred: 


Nc  < 3°  , 


Nbgi  10 ' 


N < 100,  Nr  < 10, 

(67) 

Nbr  5 10  ' Nbc<5, 

(68) 

M < 0.83  E + 75  . 

(69) 

13 


Inequality  (69)  is  characteristic  of  the  IBM  370/155  computer . 
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Here,  Nc  is  the  number  of  cycles,  H is  the  total  number  of 
iterations,  Nr  is  the  number  of  restorative  iterations  per 
cycle,  is  the  number  of  bisections  of  the  gradient  stepsize 

required  to  satisfy  Ineqs.  (59),  is  the  number  of  bisec- 

tions of  the  restoration  stepsize  required  to  satisfy  Ineqs. 

(62)  , Nj_)c  is  the  number  of  bisections  of  the  gradient  stepsize 
required  to  satisfy  Ineq.  (63),  and  M is  the  modulus  of  any  of 
the  quantities  employed  in  the  algorithm. 


7 . Numerical  Examples 

In  this  section,  eight  numerical  examples  are  described 
employing  scalar  notation.  In  particular,  the  symbols 
x^(0),i=l,...,n,  denote  the  components  of  the  original  state 
vector;  the  symbols  y^(0),i  = l,...,k,  denote  the  components  of 
the  auxiliary  state  vector;  the  symbol  u(0)  denotes  the  origi- 
nal control  variable  (a  scalar  in  the  examples  reported  here) ; 
the  symbol  w(9)  denotes  the  auxiliary  control  variable  (a  scalar) ; 
and  the  symbol  t denotes  the  final  time. 

For  all  of  the  examples,  a time  normalization  is  used  in 
order  to  simplify  the  numerical  computations.  Specifically, 
the  actual  time  0 is  replaced  by  the  normalized  time 

t = e/T  , (70) 

which  is  defined  in  such  a way  that  t = 0 at  the  initial  point 
and  t = 1 at  the  final  point.  The  actual  final  time  x,  if  it 
is  free,  is  regarded  as  a parameter  to  be  optimized.  Thus,  the 
dimension  of  the  parameter  is  p = 0 if  t is  fixed  and  p=  1 if  x 
is  free. 

The  dot  denotes  derivative  with  respect  to  the  normalized 
time , i . e . , 


= dx^/dt, 

x2  = dx2/dt, . . . , 

x = dx  /dt , 
n n ' 

(71) 

i'l  = dyj/dt, 

Y 2 ~ dy2/dt, . . . , 

yk  = dyk/dt  . 

(72) 

Concerning  the  convergence  history,  the  terminology  is  as 

follows.  N denotes  the  cycle  number,  N is  the  number  of 
c J g 

gradient  iterations  per  cycle,  is  the  number  of  restorative 
iterations  per  cycle,  N is  the  total  number  of  iterations,  P 
is  the  constraint  error,  Q is  the  error  in  the  optimality  con- 
ditions, and  I is  the  value  of  the  functional  being  optimized. 

Example  7.1.  This  example  involves  (i)  boundary  condi- 
tions of  the  fixed  endpoint  type,  (ii)  fixed  final  time  t = 1, 
and  (iii)  a linear  state  inequality  constraint  of  the  first 
order.  Problem  PI  is  as  follows: 


1=1  (x2  + u2)dt  , (73-1) 

Jo 

x^  = x2  - u , (73-2) 

x^  - 0 . 9 > 0 , (7  3-3) 

x1 (0)  = 1,  (73-4) 

x1(l)  = 1 . (73-5) 


In  this  formulation,  the  unknowns  are  the  state  variable  Xj  (t) 
and  the  control  variable  u(t). 

We  introduce  the  auxiliary  state  variable  y^t)  and  the 
auxiliary  control  variable  w(t)  defined  by 
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~ 0 • 9 = , y±  = w 


With  this  understanding,  Problem  P2  is  represented  by 


1 = 1 (x^  + u2) dt  , 


(75-1 : 


X1  = X1  “ u f 


yL  = w , 


(75-2) 


x,  - u - 2y,  w = 0 , 


(75-3) 


x2  (0)  =1  , 


y,  (0)  = /(0.1) , 


(75-4) 


xx (1)  = 1 . 


(75-5) 


In  this  formulation,  the  unknowns  are  the  state  variables  x^ (t), 
y^(t)  and  the  control  variables  u(t),  w(t). 

Since  is  linear  in  x^,  the  elimination  of  the  state 
variable  x^  is  possible  through  the  relation 


x,  = 0 . 9 + y. 


As  a consequence,  Problem  P3  is  represented  by 


1=  \ [(0.9  + y^)2  + u2]dt  , 

Jo 


(77-1) 


Y1  = w. 


(77-2) 
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(0.9+y2)2-u-  2y1w=  0 , 

(77-3) 

V 

r—i 

o 

II 

o 

r-H 

>1 

(77-4) 

Yl(l)  = /(O.l)  . 

(77-5) 

In  this  formulation,  the  unknowns  are  the  state  variable  y-^  (t) 
and  the  control  variables  u(t),  w(t).  After  Problem  P3  is 
solved,  the  function  x^(t)  can  be  computed  a posteriori  with 
(76)  . 

In  the  computer  runs,  the  following  nominal  functions  are 
employed  for  Problem  P2 : 

x1 (t)  = 1,  yx(t)  = /(0.1) , (78-1) 

u(t)  = 1,  w(t)  = 1 , (78-2) 


and  Problem  P3 : 


Yl (t)  = /(0.1) , (79-1) 

u(t)=l,  w(t)  = 1 . (79-2) 

The  numerical  results  are  given  in  Tables  1-4.  Convergence 
to  the  desired  stopping  condition  occurs  in  N=  12  iterations 
for  Problem  P2  and  N=10  iterations  for  Problem  P3.  The  CPU 
time  for  convergence  is  T=21.2  sec  for  Problem  P2  and  T=  i ; 
sec  for  Problem  P3.  The  reduction  in  CPU  time  due  to  the  use 
of  the  new  formulation  is  57 . V . 
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x1(0)=0,  x2(0)=0,  y1 (0)  = /(0.4) , (82-4) 

x±  (1)  = 1,  x 2 (1)  = 0 . (82-5) 


In  this  formulation,  the  unknowns  are  the  state  variables 
x^ (t)  , x2 (t) , (t) , the  control  variables  u(t),  w(t),  and  the 

parameter  x . 

Since  LQ  is  linear  in  x2,  the  elimination  of  the  state 
variable  x2  is  possible  through  the  relation 


x2  = 0 . 4 - 


(83) 


As  a consequence.  Problem  P3  is  represented  by 


H 

II 

H 

(84-1) 

x1  = TU  , y1  = TW  , 

(84-2) 

2 2 , n 
x^  — u — 2y^w  = 0 , 

(84-3) 

x1  (0)  = 0,  yx  (0)  = / (0.4)  , 

(84-4) 

xL(l)  = 1,  y1(l)  = ✓(0.4) . 

(84-5) 

In  this  formulation,  the  unknowns  are  the  state  variables  x^ (t), 
y^(t),  the  control  variables  u(t),  w(t),  and  the  parameter  x. 
After  Problem  P3  is  solved, the  function  x2(t)  can  be  computed 

a posteriori  with  (83). 

In  the  computer  runs,  the  following  nominal  functions  are 
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employed  for  Problem  P2: 


x1 (t)  = t , 


x2  (t)  = 0 , 


Yl(t)  = / (0.4),  (85-1) 


u (t)  = 1 , 


w(t)  = 1 , 


X = 1 , 


(85-2) 


and  Problem  P3 : 


(t)  = t , 


yx(t)  = / (0.4) , 


(86-1) 


u (t)  = 1 , 


w (t)  = 1 , 


T = 1 


(86-2) 


The  numerical  results  are  given  in  Tables  5-8.  Convergence 
to  the  desired  stopping  conditions  occurs  in  N=  34  iterations 
for  Problem  P2  and  N=  35  iterations  for  Problem  P3.  The  CPU 
time  is  T=  103.3  sec  for  Problem  P2  and  T=75.9  sec  for  Prob- 
lem P3.  The  reduction  in  CPU  time  due  to  the  use  of  the  new 
formulation  is  26.6%. 

Example  7.3.  This  example  involves  (i)  boundary  condi- 
tions of  the  fixed  endpoint  type,  (ii)  fixed  final  time  t = 1, 
and  (iii)  a linear  state  inequality  constraint  of  the  first 
order.  Problem  PI  is  as  follows: 


I = l (x2  + u2)dt  , 


;87-i ; 


x1  = x1  - u , 


(87-2) 


x^-0.8-t+t“>0  , 


(87-3) 
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xL  (0)  = 1 , 


x1 (1)  = 1 . 


(87-4) 


:37-5) 


In  this  formulation,  the  unknowns  are  the  state  variable  x^(t) 
and  the  control  variable  u(t). 

We  introduce  the  auxiliary  state  variable  y^(t)  and  the 
auxiliary  control  variable  w(t)  defined  by 


x,  - 0.8  - t + t2  = y2  , 


yx  = w 


With  this  understanding.  Problem  P2  is  represented  by 


1 = 1 (x2  + u2)  dt  , 


(89-1) 


= xx  - u , 


y i = w 


(89-2) 


X-^  — u — 1 + 2t  - 2y^w  = 0 , 


xx  (0)  = 1, 


Xj_  (1)  = 1 . 


y,  (0)  = / (0.2) , 


( 89-3 ) 


(89-4  ) 


(89-5) 


In  this  formulation,  the  unknowns  are  the  state  variables  x^ (t), 
y^(t)  and  the  control  variables  u(t),  w(t). 

Since  Lq  is  linear  in  x^ , the  elimination  of  the  state 
variable  x^  is  possible  through  the  relation 
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x-^  = 0 . 8 + t - t2  + y2  . 


As  a consequence,  Problem  P3  is  represented  by 

•1 

1=1  [(0.8  + t - t2  + y2)  2 + u2]dt  , 

0 


Yx  = w 


(90) 


(0.8  + t-  t2  + y2)2-u-  1 + 2t  - 2y1w  = 0 , 


Yl(0)  = / (0.2)  , 
Yl(l)  = 'Z  (0.2)  . 


(91-1) 

(91-2) 

(91-3) 

(91-4) 

(91-5) 


In  this  formulation,  the  unknowns  are  the  state  variable  y^ (t) 
and  the  control  variables  u(t),  w(t).  After  Problem  P3  is 
solved,  the  function  x^(t)  can  be  computed  a posteriori  with  (90). 

In  the  computer  runs,  the  following  nominal  functions  are 
employed  for  Problem  P2: 


X-j^  (t)  = 1 , 


u(t)  = 1, 


Yl(t)  = / (0.2) , 
w(t)  = 1 , 


(92-1) 

(92-2) 


and  Problem  P3: 
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The  numerical  results  are  given  in  Tables  9-12.  Conver- 
gence to  the  desired  stopping  condition  occurs  in  N = 24  iter- 
ations for  Problem  P2  and  N=  23  iterations  for  Problem  P3 . The 
CPU  time  is  T=  43.9  sec  for  Problem  P2  and  t = 29.4  sec  for 
Problem  P3.  The  reduction  in  CPU  time  due  to  the  use  of  the 
new  formulation  is  33.2%. 

Example  7.4.  This  example  involves  (i)  boundary  condi- 
tions of  the  fixed  endpoint  type,  (ii)  fixed  final  time  x = tt/2, 
and  (iii)  a linear  state  inequality  constraint  of  the  second 
order.  Problem  PI  is  as  follows: 


1=  X (U2  - 

J 0 

x^  + xt) dt  , 

(94-1) 

x1  = XU  , 

x2  = x (2  - 4x^ ) , 

(94-2) 

1 - x2  > 0 , 

(94-3) 

Xx (0)  = 0, 

x2 (0)  =0  , 

(94-4) 

x1(l)  = 1, 

x2(l)  = 0 . 

(94-5) 

In  this  formulation,  the  unknowns  are  the  state  variables  x1 (t), 
x 2 ( t ) and  the  control  variable  u(t). 

We  introduce  the  auxiliary  state  variables  y^t),  y2(t) 
and  the  auxiliary  control  variable  w(t)  defined  by 


JL ^ J 
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1 - x2  = Yl  , 


y i = Ty2  ' 


y0  = tw 


With  this  understanding,  Problem  P2  is  represented  by 


t (u2  - x2  + it) dt  , 


(96-1) 


x1  = xu,  x2  = t(2-4x1),  yi  = 'ey,,,  Y2  = TW' 


(96-2) 


4x1u  - y2  - y-jW  = 0, 


(96-3) 


X;L(0)=0,  x2(0)=0,  Yl(0)=l,  y2(0)=-l, 


(96-4) 


X;L  (1)  = 1,  x2(l)  = 0 . 


(96-5) 


In  this  formulation,  the  unknowns  are  the  state  variables  x^(t), 
x2(t),  y1('t)»  y2(t)  and  the  control  variables  u(t),  w(t). 

Since  LQ  is  linear  in  x2,  the  elimination  of  the  state 
variable  x2  is  possible  through  the  relation 

x2-l-yj.  (97) 

As  a consequence.  Problem  P3  is  represented  by 


I = \ T (U2  - X2  + Tt)dt  , 


(98-1) 


x2  = TU  , Y;l  = xy2  , 


(98-2) 


4xlU  ~ Y2  ” Y1W  = ° ' 


(98-3) 


i 
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xx (0)  = 0 , 


Yi  (0)  = 1 , 


y2(0)  = -l  , 


(98-4) 


x^l)  = 1 , 


yx(i)  = -1  . 


(98-5) 


In  this  formulation,  the  unknowns  are  the  state  variables 
x1 (t) , y^ ( t)  , y2 (t)  and  the  control  variables  u(t),  w(t). 

After  Problem  P3  is  solved,  the  function  x^(t)  can  be  computed 
a posteriori  with  (97) . 

In  the  computer  runs,  the  following  nominal  functions  are 
employed  for  Problem  P2 : 


x,  (t)  = t,  x0  ( t)  = 0 , y,  ( t)  = 1-2 1 , y0  (t)  = - 1 , 


(99-1) 


u(t)  = l/i , w(t)  = 0, 


(99-2) 


and  Problem  P3 : 


x1(t)=t,  yx  ( t)  = 1 - 2t , y2(t)=-l,  (100-1) 


u (t)  = 1/t  , w (t)  = 0 . 


(100-2) 


The  numerical  results  are  given  in  Tables  13-16.  Conver- 
gence to  the  desired  stopping  condition  occurs  in  N = 17  itera- 
tions for  Problem  P2  and  N=6  iterations  for  Problem  P3.  The 


CPU  time  is  T = 66.7  sec  for  Problem  P2  and  T = 20.0  sec  for 
Problem  P3.  The  reduction  in  CPU  time  due  to  the  use  of  the 
new  formulation  is  70.1°. 
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Example  7.5.  This  example  involves  (i)  boundary  condi- 
tions of  the  fixed  endpoint  type,  (ii)  fixed  final  time  x=l, 
and  (iii)  a linear  state  inequality  constraint  of  the  second 
order.  Problem  Pi  is  as  follows: 


I = 1 u2dt  , 

J 0 

(101-1) 

xl  = x2  ' 

x2  = u , 

(101-2) 

0 . 15  - x1  > 0 , 

(101-3) 

x1(0)  = 0, 

x2(0)  = 1, 

(101-4) 

xx(l)  = 0, 

x2 (1)  = -1 . 

(101-5) 

In  this  formulation,  the  unknowns  are  the  state  variables  x^ (t), 
X2(t)  and  the  control  variable  u(t). 

We  introduce  the  auxiliary  state  variables  y^tt),  y2(t) 
and  the  auxiliary  control  variable  w(t)  defined  by 


0.15  - xx 


(102) 


With  this  understanding,  Problem  P2  is  represented  by 

1=  \\26t, 

J 0 

X1  = X2'  X2  = U'  ^1  = ^2'  ^2  = W ' 
u + 2y2  + 2y^w  = 0 , 


(103-1) 

(103-2) 

(103-3) 
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x1(0)=0,  x2(0)=l,  y1  (0)=  /(0.15),  y2(0 ) =-l // ( 0 . 60 ),  (103-4) 

x^(l)=0,  x2(l)=-l.  (103-5) 

In  this  formulation,  the  unknowns  are  the  state  variables  x^(t), 
x2  (t)  , y1(t),  y2  (t)  and  the  control  variables  u(t),  w(t). 

Since  Lq  and  are  linear  in  x^  and  x2,  the  elimination 
of  the  state  variables  x^  and  x2  is  possible  through  the  rela- 
tions 

x1=0.15-y^,  x2  = -2y1y2  . (104) 

As  a consequence.  Problem  P3  is  represented  by 


I = \ u2dt  , 

(105-1) 

Jo 

yL  = y2  ' y2  = w ' 

(105-2) 

u + 2y2  + 2y^w  = 0 , 

(105-3) 

Y1(0)  = ✓(0-15)  , y2(0 )=-!//  (0.60)  , 

(105-4) 

yx(l)  = /(0.15)  , y2(l)=  l//(0.60)  . 

(105-5) 

In  this  formulation,  the  unknowns  are  the  state  variable  y^t), 
y2(t)  and  the  control  variables  u(t),  w(t).  After  Problem  P3 
is  solved,  the  functions  x^(t)  and  x2(t)  can  be  computed  a 
posteriori  with  (104). 

In  the  computer  runs,  the  following  nominal  functions  are 
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employed  for  Problem  P2: 

x1(t)=0,  x2(t)=l-2t,  y1(t)=/(0.15)(2t-l)2,  y2(  t)  =(  2t-l  )//(0 . 60 ),  ( 106-1 ) 

u(t)  = — (2t— l)2,  w ( t)  = 1 , (106-2) 

and  Problem  P3: 

yx(t)  = /(0.15)  (2t-  l)2,  y2  (t)  = (2t-l)//(0. 60)  , (107-1) 

u(t)  — — (2t— l)2,  w(t)=l.  (107-2) 

The  numerical  results  are  given  in  Tables  17-20.  Conver- 
gence to  the  desired  stopping  condition  occurs  in  N = 16  iter- 
ations for  Problem  P2  and  N = 17  iterations  for  Problem  P3.  The 
CPU  time  is  T=63.4  sec  for  Problem  P2  and  T=  34.7  sec  for 
Problem  P3.  The  reduction  in  CPU  time  due  to  the  use  of  the 
new  formulation  is  45.3%. 

Example  7.6.  This  example  involves  (i)  boundary  condi- 
tions of  the  free  endpoint  type,  (ii)  fixed  final  time  x=l, 
and  (iii)  a linear 
order.  Problem  PI 


state  inequality  constraint  of  the  first 
is  as  follows: 

1=  1 (x2  + x2  + u2/200)dt  , (108-1) 

Jo 

x =x2  , x2  = u-x2,  (108-2) 

8t2  -8t+1.5-x2>0, 


(108-3) 
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x1 (0)  =0, 


x2  (0)  = -1  . 


(108-4) 


In  this  formulation,  the  unknowns  are  the  state  variables  x-^(t), 
x2(t)  and  the  control  variable  u(t). 

We  introduce  the  auxiliary  state  variable  (t)  and  the 
auxiliary  control  variable  w(t)  defined  by 


8t2-8t+1.5-x2  = y2  , 


Y1  = w 


(109) 


With  this  understanding.  Problem  P2  is  represented  by 

1 


I = 


(x2  + x2  + u2/200 ) dt  , 


X1  x2 ' 


X2  = U " x2' 


yx  = w 


(110-1) 


(no-2; 


16t  - 8 - u + x2  - 2y^w  = 0 , 


(110-3) 


xx(0)=0,  x2(0)=-l,  yx(0)  = ✓(2.5).  (110-4) 


In  this  formulation,  the  unknowns  are  the  state  variables  x^(t), 
x2(t),  y-^(t)  and  the  control  variables  u(t),  w(t). 

Since  Lq  is  linear  in  x2,  the  elimination  of  the  state 
variable  x2  is  possible  through  the  relation 


x2  = 8t2  - 8t  + 1. 5 - y2  . 


(Ill) 


I 


As  a consequence,  Problem  P3  is  represented  by 


x1(0)=0,  y1  (0)  = /(2.5) . (112-4) 

In  this  formulation,  the  unknowns  are  the  state  variables  x^  { t), 
y^(t)  and  the  control  variables  u(t),  w(t).  After  Problem  P3 
is  solved,  the  function  x2(t)  can  be  computed  a posteriori  with 
(111)  • 

In  the  computer  runs,  the  following  nominal  functions  are 
employed  for  Problem  P2: 

x1(t)=0,  x2(t)=-l,  y^  (t)  = /(2 . 5)  , (113-1) 

u(t)  = 1,  w(t)  = 1 , (113-2) 

and  Problem  P3: 

x^t)  = 0,  yL(t)  = / (2.5)  , (114-1) 

u(t)  = 1,  w(t)  = 1.  (114-2) 

The  numerical  results  are  given  in  Tables  21-24.  Conver- 
gence to  the  desired  stopping  condition  occurs  in  N = 19  iter- 
ations for  Problem  P2  and  N = 19  iterations  for  Problem  P3.  The 
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CPU  time  is  T=  50.2  sec  for  Problem  P2  and  T = 36.5  sec  for 
Problem  P3.  The  reduction  in  CPU  time  due  to  the  use  of  the 
new  formulation  is  27.4%. 

Example  7.7.  This  example  involves  (i)  boundary  condi- 
tions of  the  free  endpoint  type,  (ii)  fixed  final  time  t=1, 
and  (iii)  a linear  state  inequality  constraint  of  the  second 
order.  Problem  Pi  is  as  follows: 


1=1  (x2  + x2  + u2/2O0 ) dt , 


(115-1) 


= x2  , 


*2  = u - x2  , 


(115-2) 


8t  -8t+1.5-x^>0, 


(115-3) 


x^O)  = 0 , 


x2 (0)= -1  . 


(115-4) 


In  this  formulation,  the  unknowns  are  the  state  variables 
x^ (t) , x2(t)  and  the  control  variable  u(t). 

We  introduce  the  auxiliary  state  variables  y-^(t),  y2(t) 
and  the  auxiliary  control  variable  w(t)  defined  by 


8t2  - 8t  + 1.5  - x1  = y2. 


yl  = y2  ' 


y9=  w. 


(116) 


With  this  understanding,  Problem  P2  is  represented  by 


I = \ (x2  + x2  + u2/200)  dt  , 


(117-1 : 
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*l  = x2/  x2  = u-x2,  Yi  = y 2 1 Y2  = w'  (117-2) 

16  - u + x2  - 2y2  - 2y^w  = 0 , (117-3) 

x1(0)=0,  x2(0)=-l,  y^(0)  =/ (1.5),  y2(0)  = -7//6.  (117-4) 

In  this  formulation,  the  unknowns  are  the  state  variables  x-^(t), 
x2  (t)  , y-^(t),  y 2(t)  and  the  control  variables  u(t),  w(t). 

Since  Lg  and  are  linear  in  x-^  and  x2,  the  elimination 
of  the  state  variables  x-^  and  x2  is  possible  through  the  rela- 
tions 

x^=8t2-8t+1.5-y2,  x2  = 16t-8-2y1y2.  (118) 

As  a consequence,  Problem  P3  is  represented  by 


r1 

1=  \ [(8t2-8t+l.  5-y2)  2 + (16t-8-2y  y )2  + u2/200]dt, 

(119-1) 

J 0 

yx  = y2 , y2  = w ' 

(119-2) 

16t  + 8 - u - 2y1y2  - 2y2  - 2y1w  = 0 , 

(119-3) 

Yl  (0)  = /( 1 . 5) , y2 (0 ) = -7//6  . 

(119-4) 

In  this  formulation,  the  unknowns  are  the  state  variables  y^  (t)  , 
y2(t)  and  the  control  variables  u(t),  w(t).  After  Problem  P3 
is  solved,  the  functions  x^(t)  and  x2(t)  can  be  computed  a 
posteriori  with  (118). 

In  the  computer  runs,  the  following  nominal  functions  are 
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employed  for  Problem  P2: 

x1(t)  = 0,  x2(t)  = -l,  Yl(t)=  / (1 . 5),  y2(t)=  -7//6  , (120-1) 

u ( t ) = 1,  w(t)  = 1,  (120-2) 

and  Problem  P3: 

Yl(t)  = /(1.5)  , y2(t)=-7//6,  (121-1) 

u ( t ) = 1,  w(t)  = 1 . (121-2) 

The  numerical  results  are  given  in  Tables  25-28.  Conver- 
gence to  che  desired  stopping  condition  occurs  in  N = 22  iter- 
ations for  Problem  P2  and  N = 26  iterations  for  Problem  P3.  The 
CPU  time  is  T=  81.9  sec  for  Problem  P2  and  T=  55.6  sec  for 
Problem  P3.  The  reduction  in  CPU  time  due  to  the  use  of  the 
new  formulation  is  32.2%. 

Example  7.8.  This  example  involves  (i)  boundary  condi- 
tions of  the  fixed  endpoint  type,  (ii)  fixed  final  time  t = 1, 
and  (iii)  a linear  state  inequality  constraint  of  the  third 
order.  Problem  PI  is  as  follows: 


I = \ u2dt  , 

(122-1) 

J 0 

xL  = x2  ' X2  = x3  ' 

x3  = u , 

(122-2) 

0 . 5 - Xj  >0  , 

(122-3) 
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x1(0)=0,  x2(0)=l,  x3(0)=2,  (122-4) 

x1(l)=0,  x2(l)=-l,  x3(1)=2.  (122-5) 

In  this  formulation,  the  unknowns  are  the  state  variables  x^t), 
x2(t),  x3(t)  and  the  control  variable  u(t). 

We  introduce  the  auxiliary  state  variables  (t) , y2 (t) , 
y3(t)  and  the  auxiliary  control  variable  w(t)  defined  by 


0. 5 - xx  = yx  , 


Y1  = Y2' 


Y2  = Y3  ' 


y3  = w 


(123) 


With  this  understanding,  Problem  P2  is  represented  by 

J 0 


(124-1) 


X1  = x2  ' 


x2  = x3  , 


x3  = u , 


(124-2) 


Y1  = Y2  ' 


Y2  = Y3  ' 


y3  = w 


(124-3) 


u + 6y2y3  + 2y1w  = 0 , 


(124-4) 


x1(0)=0,  x2  ( 0 ) =1  , x3(0)=2. 


(124-5) 


y2 (0)  = 1//2  , y2  ( 0 ) = -1//2  , y3 (0)  = -3//2,  (124-6) 


xL(l)  = 0 , 


x2  (1)  = -1  , 


x3(l)  = 2 . 


(124-7) 


In  this  formulation,  the  unknowns  are  the  state  variables  x^(t), 
x2(t),  x3(t),  y^(t),  y 2 ( t ) , y3(t)  and  the  control  variables 


u (t) , w(t) . 
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Since  are  linear  in  x^,  x 2,  x^,  the  elim- 

ination of  the  state  variables  x^,  x^,  is  possible  through 


the  relations 

2 

Xi=0.5-yi,  x2  = -2y1y2, 

x3  = "2y2  ' 2yly3  * 

(125) 

As  a consequence,  Problem  P3  is  represented  by 

cl 

I = \ u2dt  , 

Jo 

(126-1) 

y1  = y2  ' v2  = y3 ' y3  = w > 

(126-2) 

U + 6y2y3  + 2yLw  = 0 , 

(126-3) 

y x ( 0 ) = 1//2 , y2  (0)  = -1//2  , 

y3(0)  = -3//2, 

(126-4) 

yi(l)=l//2,  y2(l)=l//2, 

y3d)  = -3//2. 

(126-5) 

In  this  formulation,  the  unknowns  are 

the  state  variables  y 

(t). 

y2  (t) , y3 (t)  and  the  control  variables 

u (t) , w (t ) . After 

Prob- 

lem  P3  is  solved,  the  functions  x^ (t) , 

x2 (t) , x3  (t)  can  be 

computed  a posteriori  with  (125) . 

In  the  computer  runs,  the  following  nominal  functions 

are 

employed  for  Problem  P2: 

x1(t)=0,  x2  ( t ) = 1 - 2t , x3(t)  = 

= 2 , 

(127-1) 

Yl(t)=l//2,  y2(t)  = (l//2)(2t-l)  , 

y3 (t)  = -3//2, 

(127-2) 

u(t)  = 0 


w(t)  = 0 , 


(127-3) 
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and  Problem  P3: 

Yx  (t)  = 1//2 , y2(t)  = (l//2)(2t-l)  , y3(t)=-3//2,  (128-1) 

u ( t ) = 0,  w(t)  = 0 . (128-2) 

The  numerical  results  are  given  in  Tables  29-32.  Conver- 
gence to  the  desired  stopping  condition  occurs  in  N = 26  iter- 
ations for  Problem  P2  and  N = 21  iterations  for  Problem  P3.  The 
CPU  time  is  T=  171.0  sec  for  Problem  P2  and  T=60.0  sec  for 
Problem  P3.  The  reduction  in  CPU  time  due  to  the  use  of  the 
new  formulation  is  65.0%. 
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Table  1.  Convergence  history,  Example  7.1,  Problem  P2 . 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

r 

0 

0 

0 

0 

0. 14E+01 

1 

0 

3 

3 

0. 52E-09 

0. 35E+00 

1.83569 

2 

l 

2 

6 

0.15E-16 

0. 14E-01 

1.66599 

3 

l 

1 

8 

0.10E-09 

0. 24E-03 

1.65742 

4 

l 

1 

10 

0. 61E-17 

0. 15E-03 

1.65697 

5 

i 

1 

12 

0.92E-18 

0. 97E-04 

1.65678 

Table  2. 

Converged 

solution, 

Example  7.1,  Problem 

P2. 

t 

X1 

yl 

u 

w 

0.0 

1.0000 

0.3162 

1.7479 

-1.1826 

0.1 

0.9410 

0.2025 

1.3353 

-1 .1104 

0.2 

0.9095 

0.0978 

1.0097 

-0.9323 

0.3 

0.9006 

0.0247 

0.8366 

-0.5177 

0.4 

0.9000 

-0.0090 

0.8067 

-0.1866 

0.5 

0.9003 

-0.0176 

0.8104 

-0.0018 

0.6 

0.9000 

-0.0094 

0.8135 

0.1818 

0.7 

0.9005 

0.0238 

0.7863 

0.5158 

0.8 

0.9094 

0.0973 

0.6442 

0. 9396 

0.9 

0.9410 

0.2024 

0.4360 

1 . 1098 

1.0 

1.0000 

0.3162 

0.2474 

1 .1898 

= 1.00000 
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Table  3.  Convergence  history,  Example  7.1,  Problem  P3. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0. 14E+01 

1 

0 

1 

1 

0. 20E-28 

0. 76E+00 

2.00000 

2 

1 

2 

4 

0.22E-18 

0.26E-01 

1.67338 

3 

1 

1 

6 

0.10E-11 

0.24E-03 

1.65745 

4 

1 

1 

8 

0.86E-18 

0. 14E-03 

1.65697 

5 

1 

1 

10 

0.39E-19 

0.91E-04 

1.65675 

Table 

4 . Converged 

solution, 

Example  7.1,  Problem 

P3 . 

t 

X1 

Y1 

u 

w 

0.0 

1.0000 

0.3162 

1.7627 

-1 . 2060 

0.1 

0.9410 

0.2025 

1.3324 

-1.1035 

0.2 

0.9095 

0.0977 

1.0106 

-0.9379 

0.3 

0.9006 

0.0245 

0.8362 

-0.5111 

0.4 

0.9000 

-0.0081 

0.8072 

-0.1759 

0.5 

0.9002 

-0.0161 

0.8104 

0.0000 

0.6 

0.9000 

-0.0081 

0.8129 

0.1763 

0.7 

0.9006 

0.0246 

0.7858 

0.5116 

0.8 

0.9095 

0.0978 

0.6439 

0.9374 

0.9 

0.9410 

0.2025 

0.4383 

1.1037 

1.0 

1.0000 

0.3162 

0.2384 

1.2040 

T = 1.00000 
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Table  5.  Convergence  history,  Example  7.2,  Problem  P2. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0. 53E+01 

1 

0 

4 

4 

0. 69E-10 

0. 10E+00 

1.63945 

2 

l 

2 

7 

0. 70E-11 

0.15E-01 

1.59570 

3 

l 

2 

10 

0. 28E-16 

0. 74E-02 

1.59032 

4 

l 

1 

12 

0. 30E-09 

0. 51E-02 

1.58767 

5 

l 

1 

14 

0.14E-09 

0. 22E-02 

1.58626 

6 

l 

1 

16 

0.74E-11 

0. 25E-02 

1.58530 

7 

l 

1 

18 

0. 75E-11 

0. 89E-03 

1.58468 

8 

l 

1 

20 

0. 56E-12 

0. 14E-02 

1.58422 

9 

l 

1 

22 

0.71E-12 

0.40E-03 

1.58390 

10 

l 

1 

24 

0. 72E-13 

0. 87E-03 

1.58364 

11 

l 

1 

26 

0.10E-12 

0. 19E-03 

1.58346 

12 

l 

1 

28 

0. 13E-13 

0.55E-03 

1.58330 

13 

l 

1 

30 

0. 25E-13 

0.10E-03 

1.58319 

14 

l 

1 

32 

0. 19E-14 

0.33E-03 

1.58311 

15 

l 

1 

34 

0.73E-14 

0. 62E-04 

1.58304 

Table 

6.  Converged  solution. 

Example 

7.2,  Problem 

P2. 

t 

X1 

X2 

*1 

u 

w 

0.0 

0.0000 

0.0000 

0.6324 

0.8971 

-0.6363 

0.1 

0.1414 

0.1253 

0.5241 

0.8860 

-0.7298 

0.2 

0.2793 

0.2382 

0.4022 

0.8527 

-0.8070 

0.3 

0.4103 

0.3274 

0.2693 

0.7982 

-0.8704 

0.4 

0.5303 

0.3832 

0.1294 

0.7098 

-0.8595 

0.5 

0.6355 

0.3990 

0.0302 

0.6482 

-0.2694 

0.6 

0.7428 

0.3966 

0.0582 

0.6926 

0.6178 

0.7 

0.3425 

0.3595 

0.2010 

0.5483 

1.0172 

0.8 

0.9183 

0.2728 

0.3565 

0.4084 

0.9487 

0.9 

0.9712 

0.1488 

0.5011 

0.2581 

0.8747 

1.0 

1.0000 

0.0000 

0.6324 

0.1144 

0.7802 

t = 1.58304 
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Table  7.  Convergence  history,  Example  7.2,  Problem  P3. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0.48E+01 

1 

0 

5 

5 

0.41E-16 

0.14E+00 

1.65592 

2 

l 

2 

8 

0.17E-10 

0.16E-01 

1.59799 

3 

l 

2 

11 

0.11E-15 

0.97E-02 

1.59082 

4 

l 

1 

13 

0.39E-09 

0.47E-02 

1.58787 

5 

i 

1 

15 

0. 55E-10 

0. 26E-02 

1.58641 

6 

l 

1 

17 

0. 55E-11 

0. 23E-02 

1.58544 

7 

l 

1 

19 

0. 26E-11 

0 . 10E-02 

1.58481 

8 

l 

1 

21 

0. 36E-12 

0. 13E-02 

1.58434 

9 

l 

1 

23 

0. 27E-12 

0. 49E-03 

1.58401 

10 

l 

1 

25 

0. 42E-13 

0.87E-03 

1.58374 

11 

l 

1 

27 

0.44E-13 

0.25E-03 

1 . 58355 

12 

l 

1 

29 

0.73E-14 

0. 57E-03 

1.58339 

13 

l 

1 

31 

0.10E-13 

0.13E-03 

1.58327 

14 

l 

1 

33 

0.13E-14 

0.36E-03 

1 . 58317 

15 

l 

1 

35 

0. 33E-14 

0.81E-04 

1.58310 

Table 

8.  Converged  solution 

, Example 

7.2,  Problem 

P3. 

t 

X1 

X2 

yl 

u 

w 

0.0 

0.0000 

0.0000 

0.6324 

0.8974 

-0.6366 

0.1 

0.1414 

0.1253 

0.5240 

0.8362 

-0.7303 

0.2 

0.2794 

0.2383 

0.4020 

0.8530 

-0.8077 

0.3 

0.4104 

0.3276 

0.2690 

0.7980 

-0.8706 

0.4 

0.5303 

0.3831 

0.1299 

0.7078 

-0.8455 

0.5 

0.6354 

0.3988 

0.0332 

0.6490 

-0.2621 

0.6 

0.7426 

0.3963 

0.0604 

0.6920 

0.6008 

0.7 

0.8425 

0.3596 

0.2008 

0.5494 

1.0156 

0.8 

0.9183 

0.2729 

0.3564 

0.4081 

0.9492 

0.9 

0.9712 

0.1489 

0.5010 

0.2577 

0.8749 

1.0 

1.0000 

0.0000 

0.6324 

0.1158 

0.7799 

x = 1.58310 
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Table  9.  Convergence  history,  Example  7.3,  Problem  P2. 


Nc 

N 

g 

N 

r 

N 

P 

Q 

1 

0 

0 

0 

0 

0.21E+01 

1 0 

4 

4 

0.48E-14 

0.11E+01 

2.43959 

2 1 

2 

7 

0. 74E-16 

0. 90E-02 

2.10181 

3 1 

1 

9 

0. 79E-13 

0. 15E-02 

2 . 09857 

4 1 

1 

11 

0. 23E-13 

0. 25E-02 

2.09694 

5 1 

1 

13 

0. 59E-15 

0.62E-03 

2.09599 

6 1 

1 

15 

0. 99E-15 

0. 11E-02 

2.09535 

7 1 

1 

17 

0. 34E-16 

0. 30E-03 

2.09492 

8 1 

1 

19 

0. 49E-16 

0. 59E-03 

2.09460 

9 1 

1 

21 

0. 26E-17 

0. 16E-03 

2.09437 

10  1 

1 

23 

0.43E-17 

0. 34E-03 

2.09419 

11  1 

0 

24 

0. 66E-08 

0. 98E-04 

2.09396 

Table  10. 

Converged 

solution.  Example  7.3,  Problem  P2. 

t 

X1 

*1 

u 

w 

0.0 

1.0000 

0.4472 

1.3993 

-1 . 5645 

0.1 

0.9748 

0.2912 

1.0571 

-1.5565 

0.2 

0.9786 

0.1367 

0.7664 

-1.4938 

0.3 

1.0  104 

0.0222 

0.6505 

-0.6597 

0.4 

1.0400 

-0.0064 

0.8809 

-0.0527  i 

0.5 

1.0500 

-0.0070 

1.1025 

-0.0009 

0.6 

1.0400 

-0.0065 

1 .2823 

0.0540 

0.7 

1.0104 

0.0224 

1.3912 

0.6614 

0.8 

0.9786 

0.1367 

1.1495 

1.4920 

0.9 

0.9748 

0.2912 

0.8434 

1 . 5566 

1.0 

1.0000 

0.4472 

0.6007 

1.5644 

( 

; = 1 . 00000 

t 
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Table  11.  Convergence  history,  Example  7.3,  Problem  P3. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0. 16E+01 

1 

0 

3 

3 

0.18E-14 

0. 25E+01 

2.85114 

2 

l 

2 

6 

0.62E-17 

0. 18E-01 

2.10447 

3 

l 

1 

8 

0.58E-14 

0. 14E-02 

2.09848 

4 

l 

1 

10 

0.21E-15 

0.24E-02 

2.09655 

5 

l 

1 

12 

0.25E-17 

0. 52E-03 

2.09564 

6 

l 

1 

14 

0.28E-17 

0. 91E-03 

2.09511 

7 

l 

1 

16 

0. 30E-18 

0.26E-03 

2.09475 

8 

l 

1 

18 

0.17E-19 

0. 49E-03 

2.09449 

9 

l 

1 

20 

0. 32E-19 

0. 15E-03 

2.09429 

10 

l 

1 

22 

0. 30E-21 

0. 29E-03 

2.09414 

11 

l 

0 

23 

0. 63E-08 

0. 93E-04 

2.09396 

Table 

12.  Converged 

solution. 

Example  7.3, 

Problem  P3. 

t 

X1 

*1 

u 

w 

0.0 

1.0000 

0.4472 

1.3996 

-1 . 5648 

0.1 

0.9748 

0.2912 

1.0572 

-1.5571 

0.2 

0.9786 

0.1364 

0.7671 

-1 . 5001 

0.3 

1.0104 

0.0217 

0.6495 

-0.6543 

0.4 

1.0400 

-0.0063 

0.8809 

-0.0475 

0.5 

1.0500 

-0.0066 

1.1025 

-0.0011 

0.6 

1.0400 

-0.0064 

1.2823 

0.0488 

0.7 

1.0104 

0.0219 

1.3921 

0.6566 

0.8 

0.9786 

0.1365 

1.1485 

1.4981 

0.9 

0.9748 

0.2912 

0.8432 

1 . 5572 

1.0 

1.0000 

0.4472 

0.6005 

1 . 5646 

x = 1.00000 
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Table  13.  Convergence  history.  Example  7.4,  Problem  P2. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0 . 54E+01 

1 

0 

3 

3 

0.85E-08 

0.18E-01 

1.24356 

2 

l 

1 

5 

0.19E-11 

0. 55E-02 

1.23749 

3 

l 

1 

7 

0.98E-16 

0. 20E-02 

1.23593 

4 

l 

1 

9 

0.13E-16 

0.15E-02 

1.23520 

5 

l 

1 

11 

0.98E-19 

0.69E-03 

1.23481 

6 

l 

1 

13 

0.20E-18 

0.11E-02 

1.23438 

7 

l 

0 

14 

0. 14E-08 

0.24E-03 

1.23427 

8 

l 

1 

16 

0.10E-18 

0.63E-03 

1.23405 

9 

l 

0 

17 

0.27E-08 

0. 97E-04 

1.23391 

Table  14. 

Converged 

solution , 

Example 

7.4,  Problem 

P2. 

t 

X1 

X2 

yl 

y2 

u 

w 

0.0 

0.0000 

0.0000 

1.0000 

-1.0000 

0.9911 

-1 . 0000 

0.1 

0.1568 

0.3090 

0.8312 

-1.1438 

0.9914 

-0.8258 

0.2 

0.3099 

0.5875 

0.6421 

-1 . 2580 

0.9555 

-0.6201 

0.3 

0.4558 

0.8081 

0.4380 

-1.3342 

0.8903 

-0.3579 

0.4 

0.5875 

0.9494 

0.2247 

-1.3778 

0.7875 

-0.2128  | 

0.5 

0.7040 

0.9999 

0.0062 

-1.4010 

0.6968 

-0.0726 

0.6 

0.8059 

0.9541 

-0.2139 

-1.3972 

0.5969 

0.1308 

0.7 

0.8902 

0.8143 

-0.4308 

-1 .3576 

0.4717 

0.3785 

0.8 

0.9527 

0.5924 

-0.6383 

-1 . 2776 

0.3211 

0.6398 

0.9 

0.9901 

0.3108 

-0.8301 

-1.1572 

0.1517 

0.8893 

1.0 

1.0000 

0.0000 

-0.9999 

-1.0000 

-0.0262 

1.1049 

t = 1.00000 
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Table  15.  Convergence  history,  Example  7.4,  Problem  P3. 


N 


c 


N 


g 


N 


r 


N 


P 


0 


I 


0 0 0 

0 3 3 

115 
10  6 


0.79E+00 
0. 33E-09 
0.27E-16 
0.91E-10 


0.14E-02 
0. 16E-03 
0. 44E-04 


1.23414 

1.23377 

1.23373 


Table  16. 

Converged 

solution. 

Example 

7.4,  Problem 

P3 . 

t 

X1 

x2 

yl 

y2 

u 

w 

0.0 

0.0000 

0.0000 

1.0000 

-1.0000 

0.9874 

-1.0000 

0.1 

0.1555 

0.3090 

0.8312 

-1.1448 

0.9846 

-0.8395 

0.2 

0.3076 

0.5882 

0.6417 

-1.2633 

0.9542 

-0.6571 

0.3 

0.4538 

0.8098 

0.4361 

-1.3482 

0.8996 

-0.4231 

0.4 

0.5884 

0.9515 

0.2200 

-1.3969 

0.8104 

-0.1994 

0.5 

0.7075 

0.9999 

-0.0009 

-1.4115 

0.7039 

0.0126 

0.6 

0.8089 

0.9508 

-0.2216 

-1.3932 

0.5848 

0.2199 

0.7 

0.8906 

0.8090 

-0.4369 

-1.3425 

0.4537 

0.4251 

0.8 

0.9509 

0.5880 

-0.6417 

-1.2597 

0.3112 

0.6282 

0.9 

0.9879 

0.3092 

-0.8311 

-1.1453 

0.1578 

0.8279 

1.0 

1.0000 

0.0000 

-1.0000 

-0.9999 

-0.0053 

1.0214 

T = 1.57079 
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Table  17.  Convergence  history,  Example  7.5,  Problem  P2. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0. 76E+01 

1 

0 

4 

4 

0.93E-08 

0. 15E-01 

5.94696 

2 

l 

1 

6 

0. 88E-16 

0. 33E-02 

5.93350 

3 

l 

1 

8 

0.10E-15 

0.15E-02 

5.92983 

4 

l 

1 

10 

0.18E-19 

0.62E-03 

5.92838 

5 

l 

1 

12 

0.15E-18 

0.42E-03 

5.92763 

6 

l 

0 

13 

0. 40E-08 

0. 23E-03 

5.92713 

7 

l 

0 

14 

0. 34E-08 

0.13E-03 

5.92683 

8 

l 

0 

15 

0.70E-08 

0 . 16E-03 

5.92658 

9 

l 

0 

16 

0.76E-08 

0.42E-04 

5.92646 

Table  18. 

Converged 

solution , 

Example  7.5,  Problem  P2. 

t 

x2 

yl 

y 2 u w 

0.0 

0.0000 

1.0000 

0.3872 

-1.2909 

-4 .4309 

1.4169 

0.1 

0.0794 

0.6055 

0.2655 

-1.1399 

-3.4551 

1.6124 

0.2 

0.1243 

0.3091 

0.1600 

-0.9654 

-2.4738 

1.9045 

0.3 

0.1445 

0.1096 

0.0741 

-0.7394 

-1.5065 

2.7869 

0.4 

0.1497 

0.0129 

0.0162 

-0.4002 

-0.4453 

3.8483 

0.5 

0.1499 

0.0000 

-0.0040 

-0.0012 

0.0323 

4 .0484 

0.6 

0.1497 

-0.0128 

0.0160 

0.3994 

-0.4435 

3.8777 

0.7 

0.1445 

-0.1096 

0.0739 

0.7407 

-1.5101 

2.7892 

0.8 

0.1243 

-0.3092 

0.1600 

0.9659 

-2.4733 

1.8971 

0.9 

0.0794 

-0.6055 

0.2655 

1.1400 

-3.4547 

1 .6104 

1.0 

0.0000 

-1.0000 

0.3872 

1.2909 

-4 .4300 

1.4157 

U 


T = 1 . 00000 


Table  19.  Convergence  history,  Example  7.5,  Problem  P3. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0.44E+01 

1 

0 

2 

2 

0.69E-10 

0.25E+00 

6.14402 

2 

l 

1 

4 

0.10E-09 

0. 16E-01 

5.95065 

3 

l 

1 

6 

0. 28E-16 

0.41E-02 

5.93466 

4 

l 

1 

8 

0.12E-15 

0. 16E-02 

5.93041 

5 

l 

1 

10 

0.42E-19 

0.76E-03 

5.92872 

6 

l 

1 

12 

0.21E-18 

0.46E-03 

5.92786 

7 

l 

0 

13 

0.64E-08 

0. 29E-03 

5.92727 

8 

l 

0 

14 

0.44E-08 

0. 13E-03 

5.92693 

9 

l 

1 

16 

0.30E-21 

0.23E-03 

5.92676 

10 

l 

0 

17 

0. 38E-09 

0.41E-04 

5.92661 

Table  20. 

Converged 

solution. 

Example  7 . 

5,  Problem 

P3 . 

t 

X1 

X2 

Y1 

y2 

u 

w 

0.0 

0.0000 

1.0000 

0.3872 

-1.2909 

-4.4377 

1.4258 

0.1 

0.0794 

0.6054 

0.2656 

-1.1397 

-3.4540 

1.6117 

0.2 

0.1243 

0.3091 

0.1601 

-0.9655 

-2.4720 

1.8973 

0.3 

0.1445 

0.1096 

0.0740 

-0.7400 

-1.5095 

2.7963 

0.4 

0.1497 

0.0129 

0.0162 

-0.3987 

-0.4435 

3 .8649 

0.5 

0.1499 

0.0000 

-0.0038 

0.0000 

0.0308 

4 .0305 

0.6 

0.1497 

-0.0129 

0.0162 

0.3987 

-0.4435 

3.8647 

0.7 

0.1445 

-0.1096 

0.0740 

0.7400 

-1.5096 

2.7964 

0.8 

0.1243 

-0.3091 

0.1601 

0.9655 

-2.4720 

1.8973 

0.9 

0.0794 

-0.6054 

0.2656 

1.1397 

-3.4540 

1.6117 

1.0 

0.0000 

-1.0000 

0.3872 

1.2909 

-4.4377 

1.4258 

t = 1.00000 
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Table  21.  Convergence  history,  Example  7.6,  Problem  P2. 


N 

C 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0. 53E+02 

1 

0 

4 

4 

0. 97E-10 

0. 20E+00 

0.61522 

2 

l 

2 

7 

0.19E-11 

0. 22E-01 

0.33462 

3 

l 

2 

10 

0. 43E-15 

0. 46E-02 

0.23701 

4 

l 

2 

13 

0. 14E-14 

0.13E-02 

0.18428 

5 

l 

1 

15 

0.85E-10 

0.37E-03 

0.17757 

6 

l 

1 

17 

0. 20E-08 

0.41E-03 

0.17444 

7 

l 

1 

19 

0.15E-11 

0. 80E-04 

0.17283 

Table 

22.  Converged  solution 

, Example 

7.6,  Problem 

P2 . 

t 

X1 

X2 

yl 

u 

w 

0.0 

0.0000 

-1 . 0000 

1.5811 

11.8479 

-6.5926 

0.1 

-0.0529 

-0.2025 

0.9902 

4 . 0080 

-5.3566 

0.2 

-0.0603 

-0.0227 

0.4927 

-0.2990 

-4.5906 

0.3 

-0.0691 

-0.1922 

0.1105 

-2.7927 

-2.7103 

0.4 

-0.1006 

-0.4212 

-0.0347 

-2.0552 

-0.4898 

0.5 

-0.1482 

-0.5025 

-0.0509 

-0.4993 

0.0316 

0.6 

-0.1957 

-0.4203 

-0.0186 

1.2095 

0.8015 

0.7 

-0.2275 

-0.2047 

0.1572 

2.1037 

2.8354 

0.8 

-0.2384 

-0.0355 

0.5055 

0.9462 

3.7765 

0.9 

-0.2387 

0.0178 

0.8730 

0.2981 

3.5047 

1.0 

-0.2352 

0.0587 

1.2005 

0.7978 

3.0240 

r = 1.00000 
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Table  23.  Convergence  history.  Example  7.6,  Problem  P3. 


Nc 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0.70E+02 

1 

0 

4 

4 

0. 26E-10 

0. 20E+00 

0.62152 

2 

l 

2 

7 

0.79E-09 

0.12E-01 

0.29797 

3 

l 

10 

0.31E-16 

0.31E-02 

0.20454 

4 

l 

2 

13 

0. 28E-14 

0.23E-02 

0.18697 

5 

l 

1 

15 

0.13E-09 

0.32E-03 

0.17853 

6 

l 

1 

17 

0.99E-09 

0.40E-03 

0.17288 

7 

l 

1 

19 

0.33E-11 

0. 30E-04 

0.17157 

Table  24.  Converged  solution,  Example  7.6,  Problem  P3 . 


t 

X1 

X2 

yl 

u 

w 

0.0 

0.0000 

-1.0000 

1.5811 

12.8747 

-6.9174 

0.1 

-0.0503 

-0.1664 

0.9728 

3.9858 

-5.4232 

0.2 

-0.0550 

-0.0038 

0.4730 

-0.4751 

-4 . 5^48 

0.3 

-0.0628 

-0.1908 

0.1043 

-2.8834 

-2.4311 

0.4 

-0  .0943 

-0.4203 

-0.0195 

-2.0359 

-0.3992 

0.5 

-0.1417 

-0.5007 

-0.0275 

-0.4961 

0.0838 

0.6 

-0.1891 

-0.4200 

-0.0034 

1.1839 

0.5831 

0.7 

-0.2209 

-0.2022 

0.1492 

2.1949 

2.6900 

0.8 

-0.2311 

-0.0261 

0.4961 

0.9679 

3.8354 

0.9 

-0.2306 

0.0234 

0.8697 

0.2259 

3.5626 

1.0 

-0.2271 

0.0534 

1.2027 

0.6597 

3.0736 

t = 1.00000 
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Table  25.  Convergence  history,  Example  7.7,  Problem  P2. 


N 

c 

N 

g 

Nr 

N 

P 

Q 

I 

0 

0 

0 

0 

0. 37E+02 

1 

0 

4 

4 

0. 30E-15 

0.85E-01 

0.97910 

2 

l 

1 

6 

0. 19E-08 

0. 32E-02 

0.79284 

3 

l 

1 

8 

0.18E-09 

0.11E-02 

0.75471 

4 

l 

1 

10 

0.13E-14 

0 . 36E-03 

0.75132 

5 

l 

1 

12 

0. 24E-14 

0. 57E-03 

0.74892 

6 

l 

1 

14 

0.17E-15 

0.19E-03 

0.74719 

7 

i 

1 

16 

0.18E-15 

0. 34E-03 

0.74585 

8 

l 

1 

18 

0. 28E-16 

0.11E-03 

0.74482 

9 

l 

1 

20 

0. 21E-16 

0.22E-03 

0.74399 

10 

l 

1 

22 

0.55E-17 

0.77E-04 

0.74332 

Table  26. 

Converged 

solution, 

Example  7 . 

7,  Problem 

P2 . 

t 

X1 

X2 

*1 

y2 

u 

w 

0.0 

0.0000 

-1.0000 

1.2247 

-2.8577 

-3.4340 

0.8576 

0.1 

-0.1085 

-1.1435 

0.9426 

-2.7882 

-1 . 9189 

0.6509 

0.2 

-0.2256 

-1.1872 

0.6675 

-2.7060 

-1.2952 

1.0958 

0.3 

-0.3430 

-1.1382 

0.4038 

-2.5527 

0.2239 

1.9870 

0.4 

-0.4453 

-0.8539 

0.1591 

-2.3435 

3.6035 

1.7528 

0.5 

-0.5049 

-0.3254 

-0.0706 

-2.3029 

4 . 9089 

-1 . 1205 

0.6 

-0.5159 

0.0519 

-0.3098 

-2.4983 

2.1866 

-2.2304 

0.7 

-0.5044 

0.1432 

-0.5696 

-2.6831 

0.1680 

-1 . 3840 

0.8 

-0.4914 

0.1087 

-0.8434 

-2.7809 

-0.4350 

-0.6379 

0.9 

-0.4832 

0.0585 

-1.1239 

-2.8211 

-0.2957 

-0.1943 

1.0 

-0.4778 

0.0670 

-1.4063 

-2.8204 

0.8164 

0.2342 

t = 1.0000^ 
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Table  27.  Convergence  history,  Example  7.7, 


Problem  P3. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0.40E+02 

1 

0 

4 

4 

0. 49E-18 

0.89E-01 

0.98869 

2 

l 

1 

6 

0.44E-08 

0.40E-02 

0.79388 

3 

l 

1 

8 

0.11E-11 

0. 36E-02 

0.76835 

4 

l 

1 

10 

0. 24E-12 

0.81E-03 

0.75838 

5 

l 

1 

12 

0.19E-14 

0.10E-02 

0.75313 

6 

l 

1 

14 

0. 31E-14 

0. 32E-03 

0.75009 

7 

l 

1 

16 

0. 17E-15 

0.49E-03 

0.74799 

8 

l 

1 

18 

0.19E-15 

0. 17E-03 

0.74647 

9 

l 

1 

20 

0. 36E-16 

0. 29E-03 

0.74531 

10 

l 

1 

22 

0.24E-16 

0.11E-03 

0.74441 

11 

l 

1 

24 

0.62E-17 

0.18E-03 

0.74368 

12 

l 

1 

26 

0.45E-17 

0.76E-04 

0.74309 

Table  28. 

Converged 

solution , 

Example  7 . 

7,  Problem 

P3 . 

t 

X1 

X2 

yl 

y2 

u 

w 

0.0 

0.0000 

-1.0000 

1.2247 

-2.8577 

-3.4315 

0.8565 

0.1 

-0.1084 

-1.1428 

0.9426 

-2.7886 

-1.9086 

0.6433 

0.2 

-0.2254 

-1.1858 

0.6674 

-2.7075 

-1.2929 

1.0830 

0.3 

-0.3428 

-1.1381 

0.4035 

-2.5547 

0.1944 

2.0000 

0.4 

-0.4452 

-0.8562 

0.1588 

-2.3413 

3.6041 

1.8130 

0.5 

-0.5049 

-0.3244 

-0.0705 

-2.2985 

4 .9466 

-1.1524 

0.6 

-0.5158 

0.0532 

-0.3095 

-2.4985 

2.1647 

-2.2655 

0.7 

-0.5042 

0.1427 

-0.5694 

-2.6843 

0.1643 

-1.3756 

0.8 

-0.4912 

0.1093 

-0.8433 

-2.7809 

-0.4146 

-0.6262 

0.9 

-0.4828 

0.0607 

-1.1237 

-2.8205 

-0.2915 

-0.1963 

1.0 

-0.4773 

0.0654 

-1.4061 

-2.8212 

0.7183 

0.2035 

x = 1.00000 


ai  u ro  h o 


66 


AAR-137 


Table  29.  Convergence  history,  Example  7.8,  Problem  P2. 


Nc 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0 . 56E+02 

1 

0 

7 

7 

0.15E-17 

0.15E+03 

276.22402 

2 

l 

2 

10 

0.92E-21 

0 . 10E+02 

201.10661 

3 

l 

1 

12 

0.52E-11 

0.26E+01 

193.87600 

4 

l 

1 

14 

0.47E-14 

0.51E+00 

192.48949 

5 

l 

1 

16 

0.16E-16 

0.18E+00 

192.13689 

6 

l 

1 

18 

0.13E-18 

0.41E-01 

192.03952 

7 

l 

1 

20 

0. 79E-21 

0.15E-01 

192.01184 

8 

l 

0 

21 

0. 92E-08 

0. 37E-02 

192.00340 

9 

l 

1 

23 

0 . 19E-22 

0. 13E-02 

192.00105 

10 

l 

0 

24 

0.11E-09 

0.41E-03 

192.00023 

11 

l 

0 

25 

0.12E-09 

0. 10E-03 

192.00002 

12 

l 

0 

26 

0.17E-09 

0. 55E-04 

191.99993 

Table  30.  Converged 

solution.  Example  7.8,  Problem  P2. 

t 

x. 

X_ 

u 

l 

2 

3 

0.0 

0.0000 

1.0000 

2.0000 

-23.9758 

0.1 

0.1062 

1.0880 

-0.1601 

-19.2074 

0.2 

0.2112 

0.9839 

-1.8405 

-14 .3999 

0.3 

0.2981 

0.7359 

-3.0402 

-9.5955 

0.4 

0.3551 

0.3919 

-3.7596 

-4.7917 

0.5 

0.3749 

0.0000 

-3.9990 

0.0000 

0.6 

0.3551 

-0.3919 

-3.7596 

4.7918 

0.7 

0.2981 

-0.7359 

-3.0402 

9.5958 

0.8 

0.2111 

-0.9839 

-1.8404 

14 .4002 

0.9 

0.1062 

-1.0880 

-0.1601 

19.2058 

1.0 

0.0000 

-1.0000 

2.0000 

23 . 9878 

t 

y-2 

y3 

w 

0.0 

0.7071 

-0.7071 

-2.1213 

10.5894 

0.1 

0.6275 

-0.8668 

-1.0699 

10.8698 

0.2 

0.5374 

-0.9154 

0.1528 

14.1789 

0.3 

0.4492 

-0.8191 

1.8903 

21.0205 

0.4 

0.3805 

-0.5149 

4 . 2429 

23. 5208 

0.5 

0.3535 

0.0000 

5.6552 

0.0000 

0.6 

0.3805 

0 . 5149 

4 . 2428 

-23.5209 

0.7 

0.4492 

0.8190 

1.8904 

-21.0209 

0.8 

0.5374 

0.9154 

0.1528 

-14  . 1789 

0.9 

0.6275 

0.8668 

-1.0699 

-10.8685 

1.0 

0.7071 

0.7071 

-2.1213 

-10.5979 

J 


T = 1.00000 
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Table  31.  Convergence  history,  Example  7.8,  Problem  P3. 


N 

c 

N 

g 

N 

r 

N 

P 

Q 

I 

0 

0 

0 

0 

0 . 39E+02 

1 

0 

5 

5 

0. 90E-18 

0.40E+02 

220.61775 

2 

l 

1 

7 

0. 75E-08 

0.50E+01 

196 .13345 

3 

i 

1 

9 

0. 37E-15 

0. 12E+01 

193.00474 

4 

l 

1 

11 

0.15E-14 

0. 31E+00 

192.27617 

5 

l 

1 

13 

0. 10E-18 

0. 99E-01 

192.08057 

6 

l 

1 

15 

0. 45E-19 

0.27E-01 

192.02401 

7 

l 

0 

16 

0.76E-08 

0.90E-02 

192.00723 

8 

l 

1 

18 

0.94E-22 

0. 26E-02 

192.00201 

9 

l 

0 

19 

0. 60E-10 

0. 82E-03 

192.00049 

10 

l 

0 

20 

0. 23E-09 

0. 29E-03 

191.99995 

11 

l 

0 

21 

0. 24E-09 

0. 61E-04 

191.99982 

Table 

32.  Converged 

solution , 

Example  7.8, 

Problem  P3. 

t 

X1 

X2 

X3 

u 

0. 0 

0.0000 

1.0000 

2.0000 

-24.0214 

0.1 

0.1061 

1.0879 

-0.1607 

-19 .2006 

0.2 

0.2111 

0.9838 

-1.8406 

-14 . 3974 

0.3 

0.2981 

0.7353 

-3.0401 

-9.5932 

0.4 

0.3551 

0.3918 

-3.7591 

-4.7858 

0.5 

0.3749 

0.0000 

-3.9981 

0.0000 

0.6 

0.3551 

-0.3918 

-3.7590 

4 . 7858 

0.  7 

0.2981 

-0.7358 

-3.0402 

9.5932 

0.8 

0.2111 

-0.9838 

-1.8406 

14 . 3974 

0.9 

0.1061 

-1 . 0879 

-0.1607 

19.2006 

1.0 

0.0000 

-1.0000 

2.0000 

24.0213 

t 

*1 

^2 

^3 

w 

0.0 

0.7071 

-0.7071 

-2.1213 

10.6217 

0.1 

0.6275 

-0.8668 

-1.0693 

10.8672 

0.2 

0.5374 

-0.9153 

0.1533 

14 . 1786 

0.3 

0.4492 

-0.8139 

1.8907 

21.0169 

0.4 

0.3805 

-0.5148 

4.2423 

23.5029 

0.5 

0.3536 

0.0000 

5.6530 

0.0001 

0.6 

0.3805 

0.5148 

4.2422 

-23.5028 

0.7 

0.4492 

0.8189 

1 .8908 

-21.0171 

0.3 

0.5374 

0.9153 

0.1533 

-14 .1786 

0.9 

0.6275 

0.8668 

-1 . 0693 

-10.8671 

1.0 

0.7071 

0.7071 

-2.1213 

-10.6217 

t = 1 . 00000 
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8 . Discussion  and  Conclusions 

In  this  paper,  we  consider  optimal  control  problems  in- 
volving the  minimization  of  a functional  subject  to  differen- 
tial constraints,  terminal  constraints,  and  a state  inequality 
constraint.  The  state  inequality  constraint  is  of  a special 
type,  namely,  it  is  linear  in  some  or  all  of  the  components 
of  the  state  vector. 

A transformation  technique  is  introduced,  by  means  of 
which  the  inequality  constrained  problem  is  converted  into  an 
equality  constrained  problem  involving  differential  constraints, 
terminal  constraints,  and  a control  equality  constraint.  The 
transformation  technique  takes  advantage  of  the  partial  linearity  of 
the  state  inequality  constraint  so  as  to  yield  a transformed 
problem  characterized  by  a new  state  vector  of  minimal  size. 

This  concept  is  important  computationally,  in  that  the  compu- 
ter time  per  iteration  increases  with  the  square  of  the  di- 
mension of  the  state  vector. 

In  order  to  illustrate  the  advantages  of  the  new  trans- 
formation technique,  eight  numerical  examples  are  solved  by 
means  of  the  sequential  ordinary  gradient-restoration  algorithm  for 
optimal  control  problems  involving  nondifferential  constraints 
(Refs.  5-6).  Two  formulations  are  employed:  that  associated 
with  Problem  P2  and  that  associated  with  Problem  P3.  If  n 

is  the  number  of  original  state  variables  and  k is  the  order 


. - L 
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of  the  state  inequality  constraint,  the  transformed  Problem  P2 
is  characterized  by  = n + k state  variables,  while  the  trans- 
formed Problem  P3  is  characterized  by  n^  = n + b state  varia- 
bles, with  b < k. 

For  Examples  7.1  through  7.8,  the  convergence  history  and 
the  converged  solution  for  both  Problem  P2  and  Problem  P3  are 
given  in  Tables  1-32.  To  facilitate  the  comparison  between 
Problem  P2  and  Problem  P3,  summary  results  are  presented  in 
Tables  33-36. 


Table  33  shows  the  number  of  cycles  Nc , the  number  of 
gradient  iterations  EN^,  the  number  of  restorative  iterations 
ENr,  and  the  total  number  of  iterations  N for  both  Problem  P2 


and  Problem  P3. 


Table  34  shows  the  CPU  time  for  Problem  P2,  the  CPU  time 


for  Problem  P3,  the  ratio  of  CPU  times  T^/T2>  and  the  relative 
saving  in  CPU  time  (T 2'  From  Table  34,  it  appears  that 
savings  in  computer  time  ranging  from  26.6%  to  70.1?  are  obtained 
by  employing  the  formulation  associated  with  Problem  P3, 
rather  than  the  formulation  associated  with  Problem  P2. 


Table  35  shows  the  CPU  time  per  iteration  for  Problem  P2, 
the  CPU  time  per  iteration  for  Problem  P3 , the  ratio  of  CPU 
times  per  iteration  anc*  t'ie  re-'-at:'-ve  saving  in  CPU  time 

per  iteration  (x  ^ ~ T ■$)  / x2  ‘ From  Table  35,  it  appears  that 
savings  in  computer  time  per  iteration  ranging  from  15.1%  to 
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56.6%  are  obtained  by  employing  the  formulation  associated  with 
Problem  P3,  rather  than  the  formulation  associated  with  Prob- 
lem P2.  From  Tables  34-35,  the  beneficial  effects  due  to  the 
use  of  the  new  transformation  technique  are  apparent. 

Table  36  shows  the  value  of  the  functional  I at  conver- 
gence for  both  Problem  P2  and  Problem  P3.  It  is  clear  that 
the  value  of  the  functional  I obtained  with  Problem  P2  is  the 
same  as  that  obtained  with  Problem  P3  [of  course,  within  the 
tolerance  limits  specified  by  the  stopping  conditions (64 )] . 

The  converged  state  variables  and  control  variables  for 
Problem  P2  and  Problem  P3  are  also  the  same  [of  course,  within 
the  tolerance  limits  specified  by  the  stopping  conditions  (64)]. 
This  can  be  seen  by  inspection  of  Tables  1-32. 
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Table  33.  Summary  results. 


Example 

Problem  P2 

Problem  P3 

N 

c 

ZN 

g 

ZN 

r 

N 

N 

c 

ZN 

g 

ZN 

r 

N 

7.1 

5 

4 

8 

12 

5 

4 

6 

10 

7.2 

15 

14 

20 

34 

15 

14 

21 

35 

7.3 

11 

10 

14 

24 

11 

10 

13 

23 

7.4 

9 

8 

9 

17 

3 

2 

4 

6 

7.5 

9 

8 

8 

16 

10 

9 

8 

17 

7.6 

7 

6 

13 

19 

7 

6 

13 

19 

7.7 

10 

9 

13 

22 

12 

11 

15 

26 

7.8 

12 

11 

15 

26 

11 

10 

11 

21 

Table  34.  Summary  results. 


Example 

T2 

(sec) 

T 

3 

(sec) 

VT2 

<T2-T3)/T2 

7.1 

21.23 

13.32 

0.627 

0.373 

7.2 

103.34 

75.94 

0.734 

0.266 

7.3 

43.99 

29.42 

0.668 

0.332 

7.4 

66 .79 

20.03 

0.299 

0.701 

7.5 

63.44 

34.73 

0.547 

0.453 

7.6 

50.24 

36.51 

0.726 

0.274 

7.7 

81.99 

55.61 

0.678 

0.322 

7.8 

171.02 

60.00 

0.350 

0.650 
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Table  35.  Summary  results. 


Example 

T 2 

(sec/iter) 

T 

3 

(sec/iter) 

VT2 

(T2_T3)/T2 

7.1 

1.769 

1.332 

0.752 

0.248 

7.2 

3.039 

2.169 

0.713 

0.287 

7.3 

1.832 

1.279 

0.698 

0.302 

7.4 

3.928 

3.338 

0.849 

0.151 

7.5 

3.965 

2.043 

0.515 

0.485 

7.6 

2.644 

1.921 

0.726 

0.274 

7.7 

3. 726 

2.138 

0.573 

0.427 

7.8 

6.577 

2.857 

0.434 

0.566 

Table  36. 

Summary 

results 

• 

Example 

Problem  P2 

Problem  P3 

P 

Q 

I 

P 

Q 

I 

7.1 

0.92E-18 

0. 97E-04 

1.65678 

0. 39E-19 

0. 91E-04 

1 .65675 

7.2 

0.73E-14 

0. 62E-04 

1.58304 

0. 33E-14 

0.81E-04 

1.58310 

7.3 

0. 66E-08 

0. 98E-04 

2.09396 

0.63E-08 

0. 93E-04 

2.09396 

7.4 

0. 27E-08 

0. 97E-04 

1.23391 

0.91E-10 

0.44E-04 

1 .23373 

7.5 

0. 76E-08 

0 . 42E-04 

5.92646 

0. 33E-09 

0. 41E-04 

5.92661 

7.6 

0.15E-11 

0.8 OF- 04 

0.17283 

0.33E-11 

0. 30E-04 

0.17157 

7.7 

0. 55E-17 

0.77E-04 

0.74332 

0.45E-17 

0.76E-04 

0.74309 

7.8 

0.17E-09 

0.55E-04  1 

91.  99993 

0. 24E-09 

0.61E-04 

191 . 99982 
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9 . Appendix;  Computational  Effort 

At  each  iteration  of  the  gradient  phase  or  the  restora- 
tion phase  a linear,  two-point  boundary-value  problem  must 
be  solved.  If  the  method  of  particular  solutions  (Refs.  7-9) 
is  employed  in  conjunction  with  the  sequential  ordinary 
gradient-restoration  algorithm  for  optimal  control  problems 
with  nondifferential  constraints  (Refs.  5-6),  one  must  execute 
n^  + p + 1 independent  sweeps  of  the  linearized  system  governing 
the  variations  associated  with  the  gradient  phase  or  the  res- 
toration phase,  where  n^  is  the  dimension  of  the  state  vector 

14 

and  p is  the  dimension  of  the  parameter.  Note  that 

n2 = n + k for  Problem  P2,  (129-1) 

" 

n^  = b + k for  Problem  P3 . (129-2) 

Hence,  the  following  factor  is  indicative  of  the  algorithmic 
time  per  iteration: 

Wi  = 2°i  (ni  + P + 1 ) < i=2  or  3,  (130) 


14 

The  subscript  i = 2 is  employed  for  Problem  P2 , and.  t :io  .!  - 

cript  i = 3 is  employed  for  Problem  Pi.  j'l 


- 
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with  the  implication  that 

W2 = 2 (n+k) (n + k + p + 1)  for  Problem  P2,  (131-1) 

= 2 (b  + k)  (b  + k + p + 1)  for  Problem  P3  . (131-2) 

In  certain  applications , the  nondifferential  constraint 
can  be  solved  in  terms  of  the  original  control  u,  assumed  to 
be  a scalar.  Then,  one  can  arrive  at  a new  formulation  of 
Problems  P2  and  P3,  where  the  original  control  u is  absent, 
while  the  auxiliary  control  w is  present.  This  procedure  is 
feasible  with  Examples  7 . 1 , 7.3,  and  7.5  through  7.8,  owing  to 
the  fact  that  the  nondifferential  constraint  is  linear  in  u. 

It  is  not  advisable  with  Examples  7.2  and  7.4,  owing  to  the 
presence  of  square  roots  (Example  7.2)  and  singularities 
(Example  7.4)  in  the  analytical  solution  for  u. 

For  the  sake  of  argument,  assume  that  the  above  substitu- 
tion procedure  is  employed.  Then,  one  can  use  the  sequential 
ordinary  gradient-restoration  algorithm  for  optimal  control 
problems  without  nondifferential  constraints  (Refs.  13-14). 

Here,  one  must  execute  q+1  independent  sweeps  of  the  linear- 
ized system  governing  the  variations  associated  with  the  gra- 
dient phase  or  the  restoration  phase,  where  q is  the  number 
of  final  conditions.  Hence,  the  following  factor  is  indicative 
of  the  algorithmic  time  per  iteration  : 
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Wi=2ni(q+1),  i = 2 or  3,  (132) 

with  the  implication  that 

Vi 2 - 2 (n  + k)  (q  + 1 ) for  Problem  P2  , (133-1 ) 

= 2 (b  + k)  (q + 1 ) for  Problem  P3 . (133-2) 

Numerical  results  illustrating  the  above  analysis  are 
given  in  Tables  37-39.  Table  37  contains  the  principal  data 
for  Example  7.1  through  7.8,  more  specifically:  the  number  of 
original  state  variables  n,  the  order  of  the  state  inequality 
constraint  k,  the  number  of  parameters  p,  and  the  number  of 
final  conditions  q.  The  table  also  contains  the  dimension  a 
of  the  vector  , the  dimension  b of  the  vector  xD , the  trans- 
formed  number  of  state  variables  r\2  associated  with  Problem 
P2,  and  the  transformed  number  of  state  variables  n^  associa- 
ted with  Problem  P3. 

Table  38  pertains  to  SOGRA  with  nondifferential  con- 
straints (Refs.  5-6)  and  contains  the  work  factor  V>A  , i = 2 or  3, 
which  is  indicative  of  the  algorithmic  time  per  iteration  [see 
Eqs . (131)  ] . The  table  also  contains  the  ratio  W^/ w2  and 

relative  difference  (V^-  V^)/W2,  which  constitutes  an  upper 
limit  to  the  savings  in  CPU  time  per  iteration  to  be  expected 
by  using  the  new  transformation  technique  in  connection  with 
the  algorithm  of  Refs.  5-6. 
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Finally,  Table  39  pertains  to  SOGRA  without  nondifferen- 
tial constraints  (Refs.  13-14)  and  contains  the  work  factor 
, i = 2 or  3,  which  is  indicative  of  the  algorithmic  time  per 
iteration  [see  Eqs.  (133)].  The  table  also  contains  the  ratio 


W3/W2  and  the  relative  dif ference (W^- )/W^ , which  constitutes 
an  upper  limit  to  the  savings  in  CPU  time  per  iteration  to  be 
expected  by  using  the  new  transformation  technique  in  connec- 
tion with  the  algorithm  of  Refs.  13-14. 
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Table 

37. 

Principal 

data 

for  the 

examples . 

Example 

n 

k 

P 

q 

a 

b 

n2 

n3 

7.1 

1 

1 

0 

1 

1 

0 

2 

1 

7.2 

2 

1 

1 

2 

1 

1 

3 

2 

7.3 

1 

1 

0 

1 

1 

0 

2 

1 

7.4 

2 

2 

0 

2 

1 

1 

4 

3 

7.5 

2 

2 

0 

2 

2 

0 

4 

2 

7.6 

2 

1 

0 

0 

1 

1 

3 

2 

7.7 

2 

2 

0 

0 

2 

0 

4 

2 

7.8 

3 

3 

0 

3 

3 

0 

6 

3 

Table  38.  Computational  effort,  SOGRA  with  nondifferential 
constraints,  Refs.  5-6. 


Example 

W 

2 

w 

3 

VW2 

(W2-W3)/W2 

7.1 

12 

4 

0.333 

0.667 

7.2 

30 

16 

0.533 

0.467 

7.3 

12 

4 

0.333 

0.667 

7.4 

40 

24 

0.600 

0.400 

7.5 

40 

12 

0.300 

0.700 

7.6 

24 

12 

0.500 

0.500 

7.7 

40 

12 

0.300 

0.700 

7.8 

84 

24 

0.285 

0.715 

Table  39.  Computational  effort, 
constraints,  Refs.  13- 

SOGRA 

■14. 

without 

nondifferential 

Example 

w 

2 

w3 

w 

3 

/W2 

(w2-w3)/w2 

7.1 

8 

4 

0. 

500 

0.500 

7.2 

18 

12 

0. 

667 

0.333 

7.3 

8 

4 

0. 

500 

0.500 

7.4 

24 

18 

0. 

750 

0.250 

7.5 

24 

12 

0. 

500 

0.500 

7.6 

6 

4 

0. 

667 

0.333 

7.7 

8 

4 

0. 

500 

0.500 

7.8 

48 

24 

0. 

500 

0.500 
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